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Chapter  1 


Introduction 

1.1  Motivation 

A  key  aspect  in  mankind’s  exploration  of  space  has  been  the  safe  delivery 
of  human  beings  or  sensitive  equipment  to  the  surface  of  a  planetary  body.  In 
order  to  complete  this  objective,  a  vehicle  must  endure  intense  aerothermodynamic 
loads  on  its  forward  surfaces  that  must  be  absorbed  or  dissipated  using  a  thermal 
protection  system.  Since  Allen  and  Eggers  showed  the  benefits  of  using  a  blunt 
body  to  defray  the  adverse  heat  loads  of  entry  in  the  1950’s,1  all  NASA  missions 
involving  atmospheric  penetration  have  employed  some  manner  of  axisymmetric, 
blunt  heat  shield.  The  manned  capsule  missions  of  the  1960’s  and  early  1970’s 
(Mercury,  Gemini,  and  Apollo)  all  employed  a  spherical  segment,  ablating  heat 
shield  while  the  Martian  Viking  and  Pathfinder  missions  both  used  a  slightly  more 
complex  spherically  blunted  cone  geometry  to  deliver  their  hardware  to  the  surface 
of  the  Red  planet.  While  these  designs  have  served  well  in  their  various  missions,  the 
improved  ability  to  handle  heat  loads  often  comes  at  the  price  of  poor  aerodynamic 
performance,  particularly  in  regards  to  L/ D  (which  correlates  to  less  flexibility  and 
safety  in  landing).  There  is  considerable  utility  in  an  optimized  shape,  neither 
necessarily  axisymmetric  nor  entirely  blunt,  that  balances  thermodynamic  loads 
with  aerodynamic  performance  in  an  ideal  fashion. 
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To  that  end,  previous  work  at  the  University  of  Maryland  introduced  a  single 
design  point  optimization2  using  modified  Newtonian  impact  theory  to  add  newer 
shapes,  such  as  those  with  elliptical  and  polygonal  base  cross  sections,  to  the  design 
space  for  entry  vehicles.  Generally,  higher  L/D  ratios  are  generated  by  employing 
a  parallelogram  cross  section  with  more  of  a  nosecone-like  axial  profile.  Coupling 
trajectory  analysis  allowed  the  optimizer  to  explore  time  dependent  characteristics 
such  as  heat  load,  cross  range,  and  deceleration  loads  for  an  entire  entry  profile. 
Optimized  vehicle  and  trajectory  pairs  have  been  generated  using  this  method  for 
both  Lunar  return3,4  (11  km/s)  and  Mars  return4,5  (12.5  krn/s)  mission  profiles. 

The  present  work  uses  computational  fluid  dynamic  solutions  in  order  to  test 
the  underlying  assumptions  of  the  optimization  process  as  well  as  to  explore  regions 
of  the  flow  field  not  necessarily  covered  by  the  simple  surface  inclination  methods 
employed  in  their  original  derivation.  This  is  done  to  provide  a  more  robust  analysis 
of  the  conclusions  in  the  works  by  Johnson. 2-6  The  effect  of  changing  the  shoulder 
radius  and  edge  sharpness  of  optimized  geometries  on  the  resulting  aerothermo- 
dynamics  is  detailed.  Slender,  rounded  polygons  and  other  unconventional  shapes 
generated  for  earth  entry  at  extra-orbital  speeds  are  examined  in  order  to  produce 
a  clearer  picture  of  the  hypersonic  aerothermodynamic  environment  experienced  by 
these  vehicles,  and  to  determine  if  the  boundaries  of  the  current  design  space  are 
valid  or  need  to  be  altered.  Preliminary  uncoupled  radiation  solutions  are  also  ex¬ 
plored  for  select  shapes  to  understand  better  the  full  heating  profile  of  these  types  of 
geometries.  Ultimately,  the  scope  of  this  work  will  be  to  enhance  the  understanding 
of  the  design  space  for  a  lower-order  optimization  method  using  higher  fidelity  com- 
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putational  fluid  dynamics,  so  that  more  viable  results  can  be  generated  in  future 
studies. 

1.2  Previous  Work 

The  goal  of  this  work  is  to  analyze  the  merits  of  particular  heat  shield  design 
optimized  for  certain  favorable  aerothermodynamic  characteristics,  not  to  add  to 
the  design  space  for  re-entry  vehicles.  As  such,  this  section  will  not  delve  into  the 
history  of  blunt  body  entry,  but  rather  provide  a  general  overview  of  the  process 
of  Johnson,  et  al. 2-6  in  which  optimized  blunt  body  geometries  were  developed  for 
certain  design  parameters.  For  a  more  detailed  discussion  as  to  what  motivated  the 
particular  choices  that  make  up  the  optimization  model  please  see  Refs.  [5]  and 
[6].  The  aforementioned  optimization  procedure  involves  four  main  aspects:  (1) 
choosing  a  geometry,  (2)  determining  the  aerodynamics,  (3)  calculating  the  heating 
profile,  and  (4)  finding  the  optimum  balance  of  (l)-(3). 

1.2.1  Heat  Shield  Geometries 

The  heat  shield  shapes  examined  in  this  study  were  defined  by  Johnson,  et 
al.2  and  are  formed  by  sweeping  an  axial  profile  (one  of  three  different  geometric 
patterns)  around  the  axis  of  an  elliptical  base  cross  section.  The  coordinate  system 
used  in  this  work  is  shown  in  Figure  1.1.  Where  (j)  is  the  rotation  angle  of  the  cross 
section  and  u  is  a  sweep  angle  for  the  axial  profile. 
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Figure  1.1:  Fixed-body  coordinate  system  shown  on  a  60°axisymmetric  spherical 
segment  heat  shield 


1.2. 1.1  Base  Cross  Section 


The  base  is  controlled  by  Giclis’  superformula  of  the  superellipse7  with  0  < 
<  27 r  shown  as: 

n3-|  -1/ni 


r(0)  = 


'  cos(|mi0) 

n2 

+ 

sin(|mi0) 

n  3- 

V\ 

^2 

(1.1) 


This  equation  has  the  ability  to  produce  a  wide  range  of  eccentric  concave 
or  convex  with  round  or  sharp  edges  solely  by  varying  the  individual  parameters. 
Here,  the  mi  parameter  corresponds  to  the  number  of  sides  of  the  supercllipse.  All 
cases  studied  in  this  work  use  a  value  of  mi  =  4,  as  it  has  been  shown  that  this 
particular  value  produces  geometries  with  the  highest  L/D.2  For  this  value  of  mi, 
the  ri\  modifier  must  be  set  to  1  to  form  viable  designs.  In  order  to  produce  closed 
shapes,  be  they  sharp  or  round  edged,  both  iq  and  z/2  must  also  be  set  to  unity 
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and  713  must  equal  n2.  The  17,2  parameter  controls  the  concavity  and  edge  sharpness 
of  the  base.  When  n2  =  2  the  base  is  an  ellipse  (regardless  of  the  value  of  mi  or 
ni).  The  base  is  convex  when  77-2  <  2  and  concave  when  n2  >  2.  Since  convex  heat 
shields  may  be  infeasible  to  efficiently  implement,  this  work  only  considers  shapes 
with  n2  <  2.  A  sample  of  the  range  of  shapes  Eq.  1.1  is  able  to  create  by  varying 
n2  is  shown  in  Figure  1.2. 


a)  mi  =  4,  n  i  =  1.5. 


b)  it,  =  2.0. 


c)  mi  =  4,  Hi  =  4.0. 


Figure  1.2:  Range  of  shapes  produced  by  Eq.  1.1,  from  ref  [5] 


Traditional  definitions  of  eccentricity  do  not  apply  here  as  U\  —  u2  —  l.  Still, 
“eccentric”  shapes  can  be  generated  by  defining  a  new  set  of  semimajor  and  semimi¬ 
nor  axes  based  on  an  input  eccentricity  parameter,  e,  as  shown  in  the  following 
equations: 


a  2 


b2(  1-e2)^ 

1 


b2  = 


1 

< 


-1  <  e  <  0 

0  <  e  <  1 

-1  <  e  <  0 


(1.2) 


(1.3) 


a2(l  —  e2)z  0  <  e  <  1 

Here  e  is  fixed  between  -1  and  1,  and  prolate  and  oblate  shapes  are  produced  when 


e  >  0  and  e  <  0  respectively.  By  using  these  new  values  of  semimajor  and  semiminor 


5 


axes  to  scale  the  Cartesian  components  of  a  proportioned  (to  a  desired  reference 
radius)  version  of  the  superformula  of  the  superellipse,  a  wide  range  of  elliptic, 
rounded  edge  bases  may  be  produced  as  shown  in  Figure  1.3. 

0 

(a)  ?i2  =  2.0,  e  =  .75 

Figure  1.3:  Example  of  oblate  and  prolate  bases 

1.2. 1.2  Axial  Profiles 

The  heat  shield  axial  profile  was  defined  by  Johnson4,6  as  the  portion  of  the 
vehicle  that  protrudes  from  the  base.  Three  different  axial  shapes  were  used  to 
generate  the  geometries  in  the  optimizer  design  space:  (1)  a  spherical  segment,  (2) 
a  spherically  blunted  cone,  and  (3)  a  power  law.  Once  chosen,  axial  profiles  are  then 
swept  about  the  contour  of  the  base  cross  section  to  construct  the  full  3-D  geometry. 
Since  the  base  of  the  heat  shield  is  not  axisymmetric,  the  axial  profiles  need  to  be 
scaled  to  the  local  radius  of  the  base  at  a  given  particular  sweep  angle.  Thus  the 
axial  profiles  can  only  be  described  as  the  shapes  below  at  0  =  0.  More  generally, 
the  axial  shape  at  a  given  rotation  angle  is  a  scaled  version  of  the  three  classes  of 
profiles  presented  in  this  section. 
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The  spherical  segment  was  defined  as  the  section  of  sphere  encompassed  by 
a  spherical  segment  angle  6S  in  which  a  plane  parallel  to  the  yz-plane  divides  the 
sphere.  The  Apollo  Command  Module  employed  a  spherical  segment  as  its  heat 
shield,  with  0S  =  25°.  A  spherically  blunted  cone  was,  simply,  a  cone  with  its  tip 
replaced  by  a  spherical  nose.  This  profile  was  defined  by  the  cone  angle  (9C)  and 
by  the  ratio  of  nose  radius  to  base  diameter  (rn/d). An  axisymmetric  spherically 
blunted  cone  heat  shield  was  used  by  the  NASA  Viking  spacecraft  to  safely  traverse 
the  Martian  atmosphere.  Finally,  the  power  law  axial  profile  was  defined  by  the 
equation: 

y  =  Axb  (1.4) 

Where  the  A  parameter  provides  a  measure  of  bluntness  for  the  shape  while  the 
b  parameter  transforms  the  shape  from  a  sharp  cone  ( b  =  0.01)  to  a  flat  surface 
( b  =  1.0).  Figure  1.4  shows  examples  of  these  three  kinds  of  axial  profile,  ft  should 
be  noted,  however,  that  in  this  work,  no  heat  shields  with  an  axial  profile  of  a  power 
law  are  explored,  as  optimized  designs  with  that  profile  mimicked  designs  that  used 
the  other  axial  profiles. 
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(a)  Spherical  Segment  (b)  Power  Law  (.4  =  1,  b  =  0.6) 


(c)  Blunted  Cone  (rn/d  =  0.3,  9C  =  35°) 


Figure  1.4:  Axial  Shapes 


1.2.2  Aerodynamic  Model 


The  aerodynamic  characteristics  of  a  certain  design  were  calculated  based  on  a 
modified  Newtonian  surface  pressure  distribution.  Figure  1.5  shows  the  conventions 
for  a  (angle  of  attack)  and  j3  (sideslip  angle)  used  for  all  aerodynamic  calculations. 


Figure  1.5:  Freestream  coordinate  system  for  a  and  /3 


Newtonian  theory  assumes  that  component  of  a  particles  momentum  normal 
to  a  surface  is  destroyed  when  impinging  upon  that  surface,  while  its  tangential 
momentum  is  conserved. 8-10  The  pressure  coefficient,  in  Newtonian  theory,  is  given 
as: 


Cp  = 


CP.,m.ax  ( )  Vnn  ■  h  <  0 


VA  •  h  >  0 


(1.5) 


where  VA  ■  n  >  0  holds  in  vehicle’s  shadow  region,  meaning  that  the  normal  com¬ 
ponent  of  velocity  is  either  nonexistent  or  moving  away  from  the  body.  Applying 
Equation  1.5  locally,  you  get  (for  any  point  ( x,y,z )  on  the  body): 

dp  CpmaX{VXTlX  T  VyTly  T  VZTlz)  (l-®) 
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where  Vx,Vy,Vz  and  nx,ny,nz  are  the  components  of  the  free  stream  velocity  and 
local  normal  vector  respectively.  In  simple  Newtonian  theory,  it  is  assumed  that 


CPimax  =  2,  whereas  modified  Newtonian  theory  uses  the  Rayleigh  pitot  tube  for¬ 
mula11  which  relates  stagnation  pressure  after  a  shock  to  the  freestream  pressure 
as: 

po.2  n- T+yrjgA  /  u+uye  y^1  n7. 

Poo  v  7  +  1  )  V47^  -2(7-  1)/ 

This  equation,  when  normalizing  by  dynamic  pressure  (goo),  yields: 


Integrating  Equation  1.6  times  a  local  area  element  (dA)  and  a  component  of  the 
body  normal  vector (nx,ny,nz)  over  the  surface  using  Simpson’s  rule  and  then  di¬ 
viding  by  the  total  planform  area  (Ap)  provides  the  non  dimensional  normal,  axial, 
and  side  forces  ( Cn,Ca,Cy )■  These  generalized  coefficients  can  be  related  to  lift 
and  drag  by: 

CLy  =  CNcos(a)  —  CUsin(a;)  (1.9) 


Cl,h  =  Cy  cos (/3)  —  Ca  cos (a)sin(/3)  (1+0) 

Cl  =  yj  ( CL,v )2  +  ( Cl,h )2  (1+1) 

Cd  =  CArsin(a)  +  Cy  sin(/3)  +  Cacos(wv )  (1+2) 


Where  wv  is  the  wind  angle  defined  as: 


(1.13) 
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In  this  work,  analysis  is  done  on  heat  shields  at  design  points  predicted  to 
deliver  peak  instantaneous  heating.  At  these  locations  on  a  particular  trajectory, 
the  vehicle  is  traveling  at  such  high  velocities  that  it  can  reasonably  be  assumed 
that: 

vx,vy>»vz  (1.14) 

Equation  1.14  essentially  means  that,  in  general,  f3  =  0°,  or  that  there  is  no  sideslip. 
With  this  condition  and  the  x-y  plane  symmetry  of  the  geometries  (true  in  this 
work,  because  m  —  4),  the  side  force,  or  C'y  will  also  be  equal  to  0.  As  such,  lift 
and  drag  now  become: 

Cl  =  Cl,v  —  Cjycos(a)  —  Ca  sin(o:)  (1-15) 

Cd  —  Cjvsin(a:)  +  CAcos(ce)  (1-16) 

A  similar  process  was  used  to  generate  aerodynamic  moment  coefficients  for 
pitching,  rolling,  and  yawing;  but,  since  stability  analysis  is  not  performed  in  this 
work,  those  equations  are  not  included  here. 

1.2.3  Heating  Models 

The  strength  and  shape  of  a  local  bow  shock  strongly  affects  the  resulting  heat 
transfer  delivered  to  a  blunt  body  in  a  hypersonic  flow.  Since  conduction  through  a 
shock  layer  is  negligible,  only  convective  and  radiative  heat  transfer  at  the  stagnation 
point  were  considered  in  developing  the  heating  model  for  the  optimization  process. 
Convective  heat  flux  is  related  to  a  velocity  gradient  imposed  by  the  body’s  surface 
pressure  distribution,  while  radiative  heat  flux  is  controlled  by  the  thickness  of  the 
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resulting  bow  shock  layer.  The  instantaneous  heat  flux  is  defined  as  a  power  density 
in  the  form  of  heat  per  unit  area  (W / cm2)  and  can  be  integrated  along  a  trajectory,  if 
one  exists,  to  determine  a  heat  load.  Though  the  presence  of  dissociated  and  ionized 
air  in  a  hypersonic  shock  layer  will  cause  some  coupling  between  these  two  modes 
of  heat  transfer,  the  heat  transfer  model  employed  neglects  any  coupling  effects.  It 
should  also  be  noted  here  that  for  all  altitude  dependent  free  stream  quantities,  the 
1976  Standard  Atmosphere12  was  used. 

1.2.3. 1  Convection 

Convective  heat  transfer  at  a  point  is  related  to  the  gradients  of  velocity  around 
that  point,  which  are,  in  turn,  controlled  by  the  pressure  distribution.  As  shown  in 
section  1.2.2,  the  local  pressure  distribution  is  a  function  of  the  geometry  of  body 
in  the  flow.  More  specifically,  a  smaller  local  radius  of  curvature  will  generate  larger 
velocity  gradients  and,  thus,  a  greater  amount  of  heating.13 

To  account  for  stagnation  point  convective  heat  flux,  the  model  of  Tauber  and 
Menees14  was  used.  The  most  general  form  of  this  model  is: 

W  =  (1.83il0-8)r-a5(l  -  gw)p^Vl  (1.17) 

Where  gw  is  the  ratio  of  wall  enthalpy  to  total  enthalpy  (assumed  zero  here)  and  rn  is 
the  local  nose  radius  (obvious  for  a  spherical  segment  or  a  spherically  blunted  cone, 
but  some  manipulation  was  needed  to  derive  an  effective  nose  radius  for  power  law 
geometries).  This  correlation  assumes  equilibrium  flow  conditions  and  a  fully  cat¬ 
alytic  surface,  which,  in  theory,  produced  more  conservative  heat  flux  predictions.8 
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This  relationship  also  follows  the  Fay  and  Riddell15  formulation  which  purports  that 
stagnation  point  heat  flux  is  proportional  to  the  inverse  square  root  of  the  local  nose 
radius. 

1.2. 3. 2  Radiation 

Radiative  heat  flux  is  controlled  by  three  primary  factors:  (1)  nose  radius  (rn), 
(2)  shock  stand  off  distance  (the  farther  away  from  the  body  the  shock  is,  the  larger 
the  radiative  heat  transfer  will  be),13  and  (3)  angle  of  attack  (a).  These  parameters 
were  combined  to  form  an  effective  nose  radius  upon  which  to  apply  the  following 
semi-empirical  radiative  heat  transfer  relations.  For  a  sphere  at  Ro  <  9000  m/s  the 
correlation  for  radiative  heat  transfer  is: 

q.s.rad  =  reffgi(3.2 8084a;10-Voo)ff2  (1.18) 

Where  g\  =  372.6,  g2  =  8.5,  and  g3  =  1.6  for  <  7620  m/s8  and  g\  =  25.34, 
g2  =  12.5,  and  g3  =  1.78  for  7620  m/s<  V, x  <  9000  m/s.16  For  a  sphere  at 
Voo  >  9000  m/s,  the  following  relation  from  Tauber  and  Sutton1'  were  applied: 

qs,rad  =  4.376a;104rf//p^22/(Roo)  (1.19) 

where  H  =  1.072;rl06Rc^1-88p“0'325 
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where  the  value  for  f(V0 c)  was  taken  from  curve  fits  from  tabulated  values17  as: 


-3.93206793xl(T12l4  +  1.61370008a;l(rV3 


-2.4359860hrl(r3l4  +  16.107869114 


f(Vo o)  =  { 


-39494.8753 


1.0023310CtelCr12VA  +  4.8977467(M(r843 


9000m/s  <Voc<  11500  m/s 


-8.42982517a;10-4V^  +  6.25552579614 


-17168.3333 


11500ra/s  <Voo<  16000  m/s 

(1.20) 

Since  the  above  relations  are  for  spheres  and  the  optimized  geometries  can  be 
other  shapes,  the  effective  nose  radius  in  these  equations  needed  to  be  related  back 
to  the  actual  nose  radius.  To  find  this  radius,  first,  the  shock  strength  (p2/pi)  was 
calculated  using  the  method  of  Tannehill,18  which  employs  empirical  curve  fits  for 
the  specific  heat  ratio  (7)  behind  the  shock  .  Then,  the  semi-empirical  method  of 
Kaattari19,20  was  employed  to  determine  the  physical  shock  stand-off  distance  (Aso). 

In  this  method,  shock  stand  off  distance  is  related  to  the  curvature  of  the  shock  by 


the  follow  equation: 


A, 


1  +  4G  ^  -  1 


(1.21) 


Tsh  2  1  rsh 

\rn  , 

Where  G  is  determined  by  curve  fits  of  a  function  of  (p2/pi)  and  7.  The  ratio  of  rn 
to  rsh  was  found  by  manipulating  a  combination  of  further  empirical  curve  fits  and 
the  geometry  of  the  blunt  body  itself,  thus  allowing  Aso  to  be  backed  out.  Finally, 
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the  empirical  curve  fit  of  Reid21  was  applied  as  shown: 

(  (g  ~  l)2 

reff  \  pa  _  /2£2  _  i 

\  pi  V  Pi 

This  equation  defined  an  effective  radius,  which  essentially  relates  any  shape  to  an 
equivalent  sphere,  that  was  used  to  solve  Equations  1.18  and  1.19. 

In  this  thesis,  analysis  is  done  at  trajectory  points  where  peak  instantaneous 
heating  is  predicted  to  occur.  At  super-orbital  entry  velocities,  peak  heating  will 
often  occur  at  VA  >  9000m/s;  therefore,  radiation  results  in  this  study  can  only 
directly  be  related  to  values  determined  from  Equation  1.19,  the  Tauber  and  Sutton 
model. 

1.2.4  Optimization  Methods 

This  section  will  briefly  delineate  the  blunt-body  optimization  procedures  de¬ 
veloped  by  Johnson,  et  al.2-6  by  introducing  the  fundamentals  behind  the  optimiza¬ 
tion  methods,  exploring  the  objective  functions  and  constraints  used,  displaying 
sample  geometries,  and  by  discussing  the  conclusions  drawn  from  these  the  results 
for  the  two  separate  approaches  used  to  find  ideal  shapes.  For  a  robust  description 
of  the  approaches,  see  Refs.  [5]  and  [6]. 

1.2.4. 1  Single  Design  Point  Optimization 

Initially,  Johnson  sought  to  derive  optimum  blunt  body  heat  shield  geome¬ 
tries  using  a  gradient  based  algorithm!2,6  at  a  single  design  trajectory  point,  Apollo 
4  peak  heating  (h  =  61  km  and  M, x  =  32.8)  .  Optimizations  were  done  using 
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Vanderplaats  Research  &  Development,  Inc.’s  DOT  software,22  which  enacted  the 
modified  method  of  feasible  directions  (MMFD)  to  minimize  or  maximize  specific 
objective  functions  subject  to  specific  inequality  constraints  and  vehicle  design  side 
constraints.  MMFD  works  by,  first,  choosing  an  initial  feasible  design,  X1.  Then,  a 
new  (in  this  case,  second)  search  direction,  E2,  is  formulated  from  the  gradient  of 
the  objective  function  and  constraints.  A  one  dimensional  search  is  then  conducted 
to  find  a  scalar,  a*  that  will  minimize  the  particular  objective  function  in  question. 
The  scalar  is  multiplied  by  the  direction  and  added  to  the  previous  vector  of  design 
variables  to  generate  a  new  X  of  design  variables.  This  procedure  continues  until 
convergence  and  the  Kuhn- Tucker  conditions22  are  satisfied. 

Objective  functions  used  were  maximizing  Ly / D ,  maximizing  (Ly / D) / qS}tot, 
minimizing  gSjtot ,  and  minimizing  Cmfigj0l.  Optimizations  were  performed  on  each 
objective  function  for  all  three  different  choices  of  axial  profile.  The  geometric  side 
constraints  used  are  shown  in  Tabic  1.1.  Most  notably,  the  ri2  parameter  has  a 
lower  bound  of  1.3,  which  is  meant  to  prevent  the  generation  of  shapes  with  sharp 
leading  edges.  The  optimization  also  used  inequality  constraints  for  stability  and 
heat  flux  to  prevent  certain  optimizations  from  producing  entirely  infeasible  shapes. 
For  example,  when  optimizing  for  Ly / D ,  these  constraints  ensured  that  stagnation 
heat  flux  would  not  exceed  3000  W/cm2. 

Figures  1.6  and  1.7  are  examples  of  optimized  shapes  using  this  method.  Gen¬ 
erally,  this  process  showed  that  high  L/D  can  be  achieved  by  four-sided,  rounded 
edge  polynomial  cross  sections  and  that  it  was  indeed  possible  to  generate  shapes 
with  high  lift  to  drag  ratios  while  keeping  peak  stagnation  point  heating  below  1000 
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W/cm2,  less  than  what  was  required  for  NASA’s  most  recent  space  capsule  design 
(Orion).  Still,  the  design  space  used  in  this  study  often  contained  a  number  of  local 
optima,  and  optimizations  were  only  performed  at  a  single  design  point;  as  such, 
there  was  no  way  to  determine  what  vehicle  truly  represented  the  ideal. 


Table  1.1:  Design  variables  and  side  constraints  used  in  gradient  based  optimization 


Spherical  segment  Spherically  blunted  cone  Power  law 


5.0°  <9S<  89.0° 

1.3  <  n2  <  4.00 
—0.968  <  e  <  0.968 
-30°  <  a  <  30° 


55.0°  <9C<  89° 
0.15  <  rnfd  <  2.00 

1.3  <  n2  <  4.00 
—0.968  <  e  <  0.968 
-30°  <  a  <  30° 


0.900  <  A  <  10.000 

1.3  <  n2  <  4.00 
—0.968  <  e  <  0.968 
-30°  <  a  <  30° 


Figure  1.6:  Spherical  Segment  with  n2  =  1.30,  e  =  —.0968,  9S  =  89.0°,  and  a  =  18° 
optimized  for  maximum  (Lv/ D)/qSttot 


Figure  1.7:  Spherical  Segment  with  n2  =  4.00,  e  =  .0968,  9S  =  15.9°,  and  a  =  —12° 
optimized  for  minimum  qS)tot 
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1.2. 4. 2  Coupled  Vehicle/Trajectory  Optimization 

In  order  to  generate  more  robust  optimal  solutions,  a  coupled  vehicle/trajectory 
analysis  was  performed.3-5  Generating  trajectories  allowed  for  the  calculation  of 
time  dependent  parameters  such  as  heat  load  and  downrange  distance.  These  tra¬ 
jectories  were  optimized  using  the  University  of  Maryland,  College  Park  Trajectory 
Optimization  Program(UPTOP),23  which  employs  a  4th-order  Runge-Kutta  rou¬ 
tine  to  propagate  the  three-degrees-of-freedom  point  mass  equations  of  rigid-body 
motion  give  as:23-25 


dl=v 

dt 

(1.23) 

dV  B-, 

H-mFb  +  9 

(1.24) 

=  (J-^bJ)  ub  +  J-1f6 

(1.25) 

dm  (  dm  \ 

dt  \  dt  )  ^ 

(1.26) 

dqt  1  n  _ 

(1.27) 

where  p,  V,  and  g  are  given  in  an  inertial  reference  frame  and  u>,  body  rotation  rate, 
is  specified  in  a  vehicle  coordinate  system.  UPTOP  uses  a  differential  evolutionary 
scheme  (DES)26  to  find  an  ideal  solution,  where  designs,  as  in  nature,27  are  evolved 
through  generations  based  on  mutation  and  cross-over  factors  until  an  optimum  is 
found.  UPTOP  compared  favorably  to  the  Program  to  Optimize  Simulated  Tra¬ 
jectories  (POST),28  NASA’s  primary  trajectory  optimization  code,  for  an  optimal 
Space  Shuttle  ascent  trajectory  through  main  engine  cutoff.4  Multiple  objective 
functions  were  used  to  find  a  set  of  non-dominated  solutions,  or  a  Pareto  frontier, 
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to  a  given  problem.  Figure  1.8  shows  an  example  of  a  Pareto  frontier  in  which 
heat  load  ( Qs,tot )  is  minimized  and  downrange  ( Pdwn )  is  maximized  concurrently. 
Essentially,  the  Pareto  frontier  is  set  of  solutions  in  which  one  design  may  be  better 
than  another  with  respect  to  one  of  the  objective  functions,  but  not  all  of  them.27 
The  optimal  solution  lies  on  this  Pareto  frontier  and  balances  all  desired  objective 
functions  in  an  ideal  fashion 
80 
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Figure  1.8:  Pareto  frontier  example  for  a  spherical  segment  ( L/D  =  0.5,  Ve  =  12.5 
km/s) 

Optimizations  were  performed  at  both  lunar  return  {Ve  =  11  km/s)  and  Mars 
return  (Ve  =  12.5  km/s)  using  two  multi-objective  functions  sets:  (1)  maximizing 
cross  range  pxrs  (to  provide  for  more  abort  scenarios)  while  minimizing  stagnation 
point  heat  load  QS)tot  (so  that  heat  shield  mass  can  be  reduced)  and  (2)  maximizing 
down  range  pdwn  while  minimizing  stagnation  point  heat  load  Qs,tot-  Trajectories 
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were  limited  to  six  bank  angle  (0&)  changes  for  L/D  —  0.3  and  0.5  and  ten  changes 
for  L/D  =  1.0.  Trajectories  had  to  also  be  chosen  to  possess  robust  entry  corridors, 
allowable  regions  of  entry  flight  path  angle  'Je,  hi  order  to  provide  for  suitable  off- 
design  survivability  of  individual  designs.  Finally,  lower  and  upper  mass  constraints 
were  set  in  order  to  bracket  all  possible  blunt-body  solutions. 

Table  1.2  shows  the  geometric  side  constraints  used  for  this  particular  opti¬ 
mization  process.  It  should  be  noted  that  for  lunar  return,  only  spherical  segment 
geometries  were  considered  while  both  spherical  segments  and  spherically-blunted 
cones  were  examined  for  Mars  return.  The  other  axial  profiles  were  not  used  because 
optimization  using  those  shapes  produced  designs  that  mimicked  other  geometries, 
thus  not  introducing  new  heat  shields,  but  rather  reproducing  already  generated 
ones.  Another  notable  feature  of  these  geometric  constraints  is  that  the  712  sharp¬ 
ness  parameter  is  bounded  below  by  1.3,  which  allows  for  rounded  parallelogram 
cross  sections,  and  above  by  2.0,  which  prevents  any  concave  shapes  from  being 
introduced. 

Furthermore,  Table  1.3  shows  the  trajectory  and  aerodynamic  constraints 
used.  The  constraints  allowed  for  sensible  trajectories  and  for  general  aerodynamic 
stability.  Mach  number  and  altitude  final  conditions  are  based  on  suitable  values 
for  parachute  deployment,  and  deceleration  loads  ( nmax )  were  limited  to  values  less 
than  what  was  experienced  on  Apollo  10  (7^). 29,30 

Figures  1.9  and  1.10  show  examples  of  vehicles  and  corresponding  trajectories 
optimized  using  this  approach  for  both  lunar  and  Mars  return  respectively.  Gen¬ 
erally,  it  was  found  that  shapes  with  a  larger  drag  area,  specifically  CpS  (drag 
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Table  1.2:  Design  variable  constraints  used  in  vehicle/trajectory  optimization 


L/D 

Ve  =  11.0  km/s  L/D 
specific  design  variables 

Ve  =  12.5  km/s  L/D 
specific  design  variables 

0.3,  0.5 

5.0° <9S<  89.0° 

—0.968  <  e  <  0.968 
—30°  <  a<  30° 

5.0°  <6*s  <  89°,  L/D  =  0.3 

20.0°  <6*s  <  89°,  L/D  =  0.5 

55.0° <9C<  89.0° 
0.15<rn/d<  2.0 
— 0.968 <e<  0.968 
-30°  <  a <30° 

1.0 

50.0°  <9S<  89° 

— 0.968  <e<  -0.95 

0°<  a <30° 

50.0° <9S<  89° 

— 0.968<e<  -0.95 

0°<  a <30° 

Common  design  variables 

1.3  <n2<  2.00 

5s<ti  <55s 
ti  +  10s  <  t2  <  t\  +  55s 
^2  +  10s  <  ^3  <  ^2  T55 
£3  T  10s  5;  ^3  T  55s 

t4  +  10s<t5  <ti  +  3605s 
t5  +  10s  <  i4  <  ti  +  3605s 

0.27  <  L/D  <0.33 

0.47  <  L/D  <0.53 

0.95 <  L/D  <0.1.05 

0°  <  4>b,all  —  180° 

For  L/D  =  0.3&0.5 
all  =  0, 1,  2,  •  •  •  ,  5,  6 

For  L/D  =  1.0 
all  =  0, 1,2,  •••  ,9, 10 

Table  1.3:  Trajectory /aerodynamic  constraints  vehicle/trajectory  optimization 


Optimization  constraints 

Trajectory  Aerodynamic  /  Geometric 

j  —  2 

tf<  3600s 

Cm,cg,a<  -0.001 

blmax  —  5 9 

Cn,cg,f}>  0.001 

ht  <  1220/cm 

sign{CLy)Chcgfi  <  0.01 

10  km  <  htj  <  45 km 

|a|<|e  +  l°| 
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coefficient  times  reference  area),  can  decelerate  at  a  higher  altitude,  in  a  less  dense 
part  of  the  atmosphere,  and  thus  produce  the  lowest  heat  loads.  Also  increases  in 
mass  were  shown  to  correlate  as  an  almost  linear  increase  to  heat  load.  For  lunar 
return,  dramatic  improvements  in  heat  loads  and  cross  range  over  the  Orion  CEV 
design  at  L/ D  =  0.27  were  shown  by  using  5°  spherical  segment  with  a  highly  oblate 
cross  section  (e  =  —.0968)  due  mostly  to  higher  drag  area  and  L/D. 
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(a)  Geometry 


(b)  Trajectory 


Figure  1.9:  Spherical  segment  with  9S  =  75.7°,  n2  =  1.31,  and  e  =  —0.967  optimized 
for  Lunar  Return  (Ve  =  11  km/ s) 


(a)  Geometry 


Time  (s) 


(b)  Trajectory 


Figure  1.10:  Spherical  segment  with  9S  =  23.7°,  n2  =  1.66,  and  e  =  0.621  optimized 


for  Lunar  Return  (Ve  =  12.5  km/s) 
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1.3  Objectives  and  Contributions 


The  aim  of  this  thesis  is  to  analyze  the  merit  of  the  aforementioned  models 
to  truly  attain  an  optimum  blunt  body  heat  shield  solution.  A  high  fidelity  CFD 
package,  DPLR,  is  employed  to  explore  these  optimized  heat  shield  designs  in  full 
detail.  The  ability  of  stagnation  point  relations  and  semi-empirical  correlations  to 
fully  predict  the  volatile  environment  of  planetary  re-entry  are  probed  in  order  to 
see  how  well  they  stand  up  when  the  full  flow  physics  is  considered.  Ultimately,  this 
work  aspires  to  find  areas  of  possible  improvement  of  the  above  optimization  model, 
weather  it  be  in  the  constraints  and  equations  themselves,  or  by  illuminating  areas 
of  the  flow-field  not  covered  by  the  lower-order  approach  (i.e.  off  stagnation  point 
heating).  As  CFD  can  be  especially  costly,  especially  in  high  temperature  environ¬ 
ment  like  re-entry,  certain  parts  of  the  lower-order  model,  such  as  the  trajectory 
integrated  variables  of  heat  load  and  downrange,  could  not  be  analyzed  with  the 
higher  order  simulation.  Still,  the  results  presented  here  provide  invaluable  detail 
to  the  design  space  of  a  lower-order  heat  shield  optimization  process. 

Some  important  results  of  this  thesis  include:  (1)  that  the  low-order  optimiza¬ 
tion  approach  does  produce  reasonable  initial  results  for  blunt-body  aerothermo- 
dynamics  (especially  in  regards  to  the  aerodynamic  coefficients),  (2)  that  high  off 
stagnation  point  convective  heating  in  parallelogram  base  designs  implies  the  need 
to  further  constrain  the  n2  parameter  to  something  higher  than  1.3,  (3)  that  the 
semi-empirical  relation  used  to  predict  convective  and  radiative  heat  transfer  break 
down  when  used  on  more  eccentric  shapes  rather  than  just  spheres,  and  (4)  that, 


24 


when  practical  matters  are  considered,  a  simple  axisymmetrical  spherical  segment 
blunt  body  heat  shield  provides  near  optimal  performance  for  both  lunar  and  Mars 
return  without  the  need  to  manufacture  exotic  shapes.  These  results  come  out  of 
extensive  benchmarking  of  the  computational  tools  (Chapter  2)  using  the  Apollo 
capsule  at  Apollo  4  peak  heating  conditions  as  a  test  subject  (Chapter  3).  A  ge¬ 
ometry  from  the  single  point  optimization  procedure2,6  is  then  examined  (Chapter 
4)  parametrically  to  fully  understand  the  effects  of  changing  certain  geometric  pa¬ 
rameters.  Chapter  5  looks  at  optimized  shapes  for  lunar3,4  and  Mars4,5  return,  and 
explore  how  these  designs  hold  up  in  a  fully  realized  flow.  Finally,  Chapter  6  con¬ 
tains  gathers  all  important  conclusions  in  full  detail  and  provides  a  summary  of  this 
work’s  important  contributions  to  the  state  of  the  art. 
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Chapter  2 


Computational  Tools 

This  chapter  provides  an  overview  of  the  various  computational  resources  em¬ 
ployed  to  simulate  blunt  body  re-entry  flows.  Only  pre-existing  software  was  used; 
however,  modifications  were  necessary,  especially  to  ensure  that  the  flow  solver  would 
properly  function  on  the  specific  hardware  found  at  the  University  of  Maryland, 
College  Park.  Essentially,  a  blunt-body  CFD  run  involves  taking  the  geometric  pa¬ 
rameters  from  the  optimizer  (as  described  in  Chapter  1),  creating  a  volume  mesh 
from  those  parameters,  converging  a  hypersonic  solution  using  the  CFD  flow  solver, 
post-processing  to  examine  the  results,  and  then,  finally,  calculating  the  thermal 
effects  due  to  shock  layer  radiation  (if  necessary).  This  general  work-flow  pattern  is 
shown  in  Figure  2.1.  The  dotted  line  connecting  DPLR,  the  flow  solver,  to  NEQA1R, 
the  radiation  solver,  signifies  that,  though  radiation  and  convection  are  considered 
uncoupled  in  this  work,  there  is  a  process  by  which  the  two  heating  regimes  can 
be  loosely  coupled.  The  individual  components  used  in  this  work  (in  the  rectangles 
outside  the  dotted  box)  are  described  in  the  following  sections. 
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Figure  2.1:  Generalized  work-flow  pattern  for  blunt-body  simulations 


2.1  Flow  Solver 


The  high  temperature  environment  of  re-entry  is  modeled  using  the  Data  Par¬ 
allel  Line  Relaxation  (DPLR)  code31  developed  at  NASA  Ames  Research  Center. 
DPLR  is  a  Fortran  90  structured  multiblock  hypersonic  continuum  code  that  utilizes 
the  Message  Passage  Interface  (MPI)  to  spread  its  workload  over  multiple  computer 
processors.  DPLR  was  chosen  because  of  its  ability  to  accurately  and  efficiently 
model  blunt  body  re-entry  at  orbital32  and  extra-orbital  velocities,33  and  because  it 
was  far  better  suited  to  the  hardware  found  at  the  University  of  Maryland,  College 
Park  than  other  hypersonic  continuum  solvers,  like  LAURA.34,35  Both  DPLR  V3.05 
and  V4.0  are  used  in  this  work.  There  is  backwards  compatibility  between  the  two 
versions  as  no  changes  were  made  to  the  parts  of  the  flow  solver  used  in  this  work 
during  the  upgrade. 

2.1.1  Nonequilibrium  Flow  Model 

The  flow  is  modeled  in  DPLR  using  the  chemically  reacting  conservation  equa¬ 
tions  equations  (derived  assuming  continuum  flow  and  that  translational  tempera¬ 
ture  T,  vibrational  temperature  Ty,  and  electron  temperature  Te  are  all  different) 
shown  in  general  form  as:36,37 
a)  Species  Continuity 

=  ^pspU° +  f>spV*p)  =  'dJsp  (2-1) 

where  sp  is  an  individual  species,  xJ  is  jth  component  of  the  orthogonal 
coordinate  directions,  wsp  is  the  mass  production  rate  of  a  species  sp  per 


unit  volume,  psp  is  the  mass  density  of  species  sp,  u30  is  the  jth  component 
of  mass  averaged  velocity,  and  Vj  is  the  jth  component  of  the  species 
mean  velocity. 


b)  Overall  Continuity 


dp  d  . 

Ft  +  aF(pu3)  =  0 


(2.2) 


where  p  is  the  sum  of  the  individual  species  densities. 


c)  ith  Direction  Species  Momentum  Conservation 


d  d  d 

~r^(PspUsp)  +  ~Q^j(PspU°Uo)  'EM  ^PSPU°^SP  d"  PspUiVsp) 


dt 


dx3 


d Psp  dr% 


dxi 


_  __  rpi  I  rpi 

Qy>j  r  elastic, sp  '  r  electric,  sp 


(2.3) 


where  t1'3  is  the  shear  stress,  Feiectric  is  the  electric  force  (a  function 
of  electric  held),  and  Feiastic  is  the  force  generated  by  elastic  collisions 
between  molecules  (from  kinetic  theory). 

d)  ith  Direction  Overall  Momentum  Conservation 

d  d  dp  dr 13  v — 

Ft  +  a? MK)  +  £>~FF  =  F.  CW  (2-4) 

sp 

where  r1’3 ,  p,  and  Feiastic  are  now  summed  over  all  species. 

e)  Electron  (denoted  by  e  subscript)  Energy  Conservation 


d_ 

dt 


Pe  fU20  +  Ccyi  +  e 


+ 


_d_ 
dx 3 


PeK  (  ^u20  +  ee 


+  Ff  +  L  +ft”Xt)  + 


d 


dx3 


id^cFe  ^  )  Feiectric, sp  Q elastic,  sp  “h  Q inelastic, 


sp 


(2.5) 
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where  ee  is  the  energy  per  unit  mass  of  an  electron,  Peiectric  is  a  power 
supplied  by  an  electric  held  on  charged  particles,  Qeiastic  is  the  rate  of 
energy  due  to  elastic  collisions,  and  Q inelastic  is  the  energy  change  due  to 
inelastic  collisions  (including  radiation). 


f)  Species  Vibrational  Energy  Conservation 


d_ 

dt 


{pSjVVyS'p) 


C) 

dxi 


(PspCV, 


v? N 

Sp^o; 


d  (  ,  dTy\ 
dxP  \nv’sp  dxp  ) 


() 

dxi 


{PspeV,spVSp) 


TT-V,sp 


+  Psp 
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where  ey:Sp  is  the  vibrational  energy  for  a  species  sp,  n'Vsp  is  the  thermal 
conductivity  for  vibrational  energy,  tt-viSP  and  re_y)Sp  are  the  T-V  and 
e-V  relaxation  times  for  a  species  sp  respectively,  and  e*V,sP  and  e*y*sP  are 
the  equilibrium  vibrational  energy  of  species  sp  at  T  and  Te  respectively. 


g)  Overall  Energy  Conservation 
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where  e  is  the  total  thermodynamic  energy  per  unit  mass,  qi  is  the  j 


th 


component  of  the  overall  heat-flux  vector,  and  Qrad  is  the  radiation  loss. 


The  above  equations  are  further  simplified,  by  choosing  the  appropriate  sim¬ 
plifications  within  DPLR  itself.  Due  to  the  high  velocities  experienced  upon  Lunar 
and  Mars  return,  an  11  species  (N2,  02,  NO,  NO+,  02+,  N,  O,  N+,  0+,  e),  19 
reaction  finite  rate  chemistry  model  for  air  from  Park38  is  considered  in  order  to 
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capture  the  effects  of  ionization  caused  by  thermal  and  chemical  nonequilibrium 
expected  after  the  bow  shock.  The  flow  is  assumed  to  be  in  thermal  nonequilibrium 
using  the  two-temperature  model  of  Park39  which  uses  an  averaged  temperature, 
defined  as: 

Td  =  VT%  (2.8) 

to  control  dissociation  rates.  Ionization  reactions  are  governed  by  bulk  translational 
temperature,  T,  as  in  the  work  by  Olynick  et  al.,40  and  translational  and  vibrational 
energy  modes  are  modeled  by  a  Landau- Teller  formulation,  which  uses  relaxation 
times  from  Milikan  and  White.41  Viscous  transport  and  thermal  conductivity  are 
modeled  using  the  mixing  rules  of  Gupta  et  ah, 42  while  species  diffusion  coefficients 
are  calculated  using  the  self-consistent  effective  binary  diffusion  (SCEBD)  method 
of  Ramshaw.43  Only  the  three  dimensional  laminar  versions  of  the  governing  equa¬ 
tions  are  considered  in  this  work,  which  is  a  reasonable  assumption  for  blunt-bodies 
entering  Earth’s  atmosphere.  Since  many  different  materials  exist  for  use  in  thermal 
protection,  a  super-catalytic  boundary  condition  is  used  for  the  heat  shield  surface. 
A  super-catalytic  surface  assumes  that  the  chemical  composition  of  the  body  is  iden¬ 
tical  to  that  in  the  freestream,  resulting  in  conservative  heating  estimates  useful  for 
design  studies.  Consequently,  this  also  means  that  material  response  is  neglected, 
as  no  specific  surface  material  is  chosen.  The  heat-shield  surface  is  also  assumed  to 
be  in  radiative  equilibrium,  in  which  energy  incident  to  the  surface  is  radiated  back 
into  the  freestream  based  on  the  equation: 

qw  =  ecrT *  (2.9) 
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where  e  is  the  surface  emissivity  (a  constant  0.85  in  this  work)  and  o  is  the  Stefan- 
Boltzmann  constant.  This  model  provides  accurate  heating  predictions,  especially 
for  the  non-ablating  heat-shields  explored  in  this  work.  Various  other  boundary  con¬ 
ditions,  such  as  periodic  or  symmetric,  can  be  employed  at  other  mesh  boundaries, 
but  the  flow  at  surfaces  through  which  air  exits  must  be  supersonic.  Due  to  limi¬ 
tations  in  the  DPLR  software,  shock  layer  radiation  is  neglected  in  the  flow  solver; 
however,  uncoupled  contributions  to  the  total  surface  heat  flux  from  the  radiating 
shock  layer  are  calculated,  in  some  cases,  with  the  NEQA1R  software  package  (see 
Section  2.3). 

2.1.2  Numerical  Model 

DPLR  is  a  fully  three-dimensional  implicit,  upwind  Navier  Stokes  solver  that 
takes  into  account  the  physical  models  discussed  above.  Euler  fluxes  are  computed 
using  a  modified  form  of  Steger- Warming  flux  vector  splitting  method  developed  by 
MacCormack  and  Candler,44  which  has  less  dissipation  than  the  original  scheme. 
Third  order  spatial  accuracy  is  maintained  by  a  MUSCL  (Monotone  Upstream- 
centered  Schemes  for  Conservation  Laws)  extrapolation  with  a  minmod  limiter.45 
A  central  differencing  approach  is  used  to  ensure  second  order  accuracy  of  the 
viscous  fluxes.  Time-marching  is  achieved  with  the  data-parallel  line  relaxation 
method(DPLR),31  which  gives  the  flow  solver  its  name.  The  DPLR  method  is  es¬ 
sentially  a  modified  version  of  McCormack’s  Gauss-Seidcl  line  rclaxation(GSLR)46 
in  that  it  uses  line  relaxation  steps  instead  of  Gauss-Seidel  sweeps  for  more  efficient 
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parallelization. 

The  best  way  to  illustrate  how  the  DPLR  method  works  is  by  applying  it  to 
the  two-dimensional  fully  implicit  inviscid  form  of  the  Navier  Stokes  equations,31 
defined  as: 


Un+ 1  -  Un  d Fn+1  d Gn+1 


=  0 


(2.10) 


At  <9£  dr] 

where  n  is  a  time  step,  U  is  the  vector  of  conserved  quantities,  and  F  and  G  are  the 

flux  vectors  in  the  body  normal  (£)  and  body  tangential  (rj)  directions.  The  fluxes 
can  be  linearized  by: 
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where  A  and  B  are  the  Jacobian  matrices  of  the  flux  vectors.  The  fluxes  are  now 
split  based  on  the  sign  of  the  eigenvalues  of  the  Jacobian  matrix  as: 


F+  +  =  A+U  +  A_U  =  F 

G+  +  G-  =  B+U  +  B_U  —  G  (2.12) 

which  allows  Equation  2.10  to  written  in  an  upwind  finite  volume  form  as: 


m,  +  (Ai/Vy)[(yl+i+jJS(+./£7y  - 
-(A+yS^jSUu  -  A_i+iJSi+ySU,-,j) 

HB_i_ySl_ij6Uij-B_i+ySi+y6Ui-1J)  = 


(2.13) 


where  is  the  solution  change  due  to  fluxes  at  time  n,  Sij  is  the  surface  area  of 
face  i,j,  and  VtJ  is  the  cell  volume.  The  DPLR  method  is  formed  by,  first,  moving 
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the  body  normal  terms  in  Equation  2.13  to  one  side,  resulting  in: 


BijSUij+i  +  Aj  j5U, 
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where  the  hatted  matrices  are  defined  as: 


4j  =  /  +  (A«/U,)(^+j+.,Si+., 


A-i-ySt-tj  +  B+i+ijSl+,j 


Bij  =  (^tlV,j)(B_,+yS,+y) 

Ci,i  =  (Ai/Vry)(B+j_ijSi_iij) 

Sy  =  (At/KyMA^.^^)  (2.15) 


Then,  the  /cmaa.  line  relaxation  steps  are  applied,  by  first  solving  a  block  tridiagonal 
system,  formed  by  neglecting  the  implicit  terms  in  2.14,  for  511^  as: 


BijSU^  +  S, 
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then,  for  /.'=!: 
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8U?j  =  Ufj™  (2.17) 

Essentially,  this  method  requires  one  LU  substitution  (for  solving  Equation  2.16) 
and  krnax  +  1  back  substitutions  for  a  single  iteration  in  time.  The  hatted  matrices 
only  need  to  be  calculated  once  for  each  relaxation  sweep,  and  each  relaxation  step 
can  be  done  in  parallel  if  body  normal  information  is  stored  locally. 
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In  the  DPLR  software  itself,  the  above  method  is  applied  to  the  fully  viscous 
three-dimensional  Navier-Stokes  equations,  which  makes  the  equations  slightly  more 
complicated;  but,  the  fundamentals  behind  the  method  are  the  same.  Convergence, 
in  DPLR,  comes  when  the  L2  norm  of  a  conserved  variable  (in  this  work,  total  density 
p  is  used)  reaches  a  sufficiently  low  level,  which  essentially  means  the  solution  is  not 
changing  between  time  steps  and  has  reached  a  quasi-steady  state.  DPLR  has  the 
ability  to  simulate  an  unsteady  flow,  but  since  the  heat  shields  studied  in  this  work 
are  examined  at  specific  trajectory  points,  this  feature  is  not  used. 

2.1.3  Work- flow 

A  typical  DPLR  simulation  for  a  given  blunt-body  heat  shield  case  is  given  as 
follows.  First,  a  Plot3D  volume  grid  is  read  into  the  software  using  a  built  in  pre¬ 
processing  package,  called  f convert .  This  package  converts  the  mesh  into  something 
that  DPLR  can  understand  while  also  breaking  it  down  into  smaller  blocks.  The 
size  and  number  of  these  blocks  are  chosen  by  the  user  to  divide  the  computation 
as  evenly  and  as  efficiently  as  possible  over  the  available  computational  nodes  to 
ensure  rapid  convergence.  Essentially,  whichever  node  has  the  most  amount  of  work 
assigned  to  it  will  dictate  the  convergence  time  of  a  given  run.  Once  the  grid  is 
sufficiently  divided,  DPLR  itself  is  invoked,  using  options  supplied  in  an  input  hie, 
which  contains  a  CFL  number  schedule  to  control  the  time  steps  of  the  simulation. 
A  first  run  of  DPLR  will  often  be  run  on  an  un-adapted  mesh  using  a  simpler 
gas  model  (i.e.  one  with  5  species  instead  of  11)  in  order  to  quickly  generate  a 
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baseline  solution.  Once  this  initial  solution  is  converged,  the  mesh  is  smoothed 
using  the  built-in  adaption  techniques  of  Saunders,47  which  reshapes  a  volume  grid 
based  on  Mach  number  contours.  Then,  the  more  complicated  gas  model  is  applied, 
and  DPLR  is  run  a  few  more  times  (usually  four  or  five)  to  convergence,  adapting 
the  grid  before  each  run.  Once  the  solution  has  reached  a  quasi-steady  state,  it 
is  considered  converged  (see  previous  section).  Finally,  a  built-in  post-processor, 
called  Postflow  is  used  to  extract  the  pertinent  flow  characteristics  on  any  part  of 
the  volume  mesh  for  view  in  any  graphical  interpreter,  such  as  Tecplot.  Postflow  also 
has  the  ability  to  integrate  flow  variables  over  surfaces,  which  is  useful  in  generating 
non-dimensional  aerodynamic  coefficients  for  lift  and  drag  (as  long  as  the  proper 
reference  quantities  are  defined). 

2.1.4  Code  Modifications  and  Validation 

DPLR  was  written  primarily  for  use  on  the  Columbia  supercomputing  cluster 
at  NASA  Ames  Research  Center;  as  such,  it  was  necessary  to  alter  the  software’s 
source  code  for  proper  function  on  a  cluster  at  LIniversity  of  Maryland,  College 
Park  (more  detail  on  both  clusters  can  be  found  in  Section  2.4).  The  changes,  while 
not  trivial,  amount  mostly  to  differences  in  semantics  used  by  various  Fortran  90 
compilers.  That  being  said,  validation  cases  were  run  using  the  modified  version  of 
DPLR  (the  changes  were  identical  for  both  DPLR  version  3.05  and  DPLR  version 
4.0)  using  a  series  of  sample  hies  distributed  with  the  DPLR  software  package.  Table 
2.1  shows  the  different  compiler /architecture  sets  used  in  this  validation  study. 
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All  sample  cases  generated  the  exact  same  results  for  all  the  compilers/architectures 


used.  For  example,  figures  2.2  and  2.3  show  Mach  contours  on  the  symmetry  plane 
for  the  Mars  Science  Laboratory(MSL)  and  surface  skin  friction  coefficients  for  a 
2-D  cylinder  respectively.  Essentially,  these  results  suggest  that  no  errors  were  in¬ 
troduced  into  the  code  by  the  modifications  made  to  the  source. 


Table  2.1:  Compilers  and  architectures  used  in  DPLR  validation  study 


Compiler 

MPI  Package 

Platform 

Cluster  Name 

Pathscale48 

Open  MPI49 

AMD  64  bit 

Skystreak  (UMD) 

Intel  Fortan50 

Open  MPI 

Intel  32  bit 

Columbia 

Gnu-Fort  an51 

Open  MPI 

Intel  32  bit 

- 

Portland  Group  Fortan52 

LAM/MPI53 

AMD  64  bit 

- 
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PATHSCALE-AMD64(OPENMPI) 


GFORTRAN-IA32  (OPENMPI) 


INTEL-IA32  (OPENMPI) 


Figure  2.2:  Mach  number  contours  on  symmetry  plane  for  MSL 


Figure  2.3:  Surface  skin  friction  coefficients  on  a  2-D  cylinder 


2.2  Grid  Generation 


Topologies  are  created  using  the  commercially  available  elliptic  grid  generation 
package,  Pointwise.54  This  package  is  an  upgrade  to  the  commonly  used  Gridgen55 
software  in  the  sense  that  it  has  a  more  streamlined  graphic  user  interface  (GUI) 
and  enhanced  undo  capabilities.  The  later  feature  is  especially  useful,  since  grid 
generation  will  often  devolve  into  a  trial-and-error  procedure  where  multiple  meshes 
are  created  until  one  “looks”  right  A  good  “looking”  mesh  will  often  have  orthogo¬ 
nality  near  its  the  boundaries  and  have  no  adverse  stretching  in  its  cells).  Pointwise 
can  create  both  structured  and  unstructured  meshes,  but  only  structured  grids  can 
be  input  to  DPLR.  Once  a  surface  grid  is  created,  Pointwise  can  improve  the  quality 
of  the  mesh  by  using  different  control  functions  (Laplace,  Middlecoff-Thomas,56  and 
Steger-Sorenson57)  to  iteratively  solve  Poisson’s  elliptical  partial  differential  equa¬ 
tions  given,  in  the  computational  domain,  by:58 

ax  ft  -  2  /3x£V  +  7W  =  -J2(Px£  +  Qxv) 

ay#  -  2 Pyft  +  7 m  =  ~j2(py^  +  QVv)  (2-18) 

with: 

a  =  x2  +  y2  f}  =  x^xv  +  y^yv  7  =  x\  +  y\ 

where  (x,  y)  are  the  Cartesian  coordinates  of  the  mesh,  (£,  rj)  are  the  transformed 
mesh  points  in  the  computational  domain,  J  is  the  Jacobian  of  the  transformation 
from  Cartesian  to  computational  domain  (J  =  x^y^  —  x^y^),  and  P  and  Q  are  source 
terms  that  provide  control  over  internal  mesh  spacing.  Pointwise  also  employs  vari- 
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ous  methods  of  hyperbolic  extrusion  in  order  to  generate  smooth  three  dimensional 
volume  meshes.  Simple  normal  extrusion  is  often  sufficient  for  generating  three- 
dimensional  blnnt-body  topologies. 

2.2.1  Work- flow 

The  process  for  creating  a  blnnt-body  mesh  for  an  optimized  heat  shield  is 
as  follows.  First,  a  Matlab  script  is  used  to  generate  the  two  profiles  (axial  and 
base  cross-section)  defined  in  Section  1.2.1  based  on  the  geometric  parameters  given 
by  the  optimizer.  These  profiles  are  fed  into  a  computer  aided  drafting  (CAD) 
software  package,  called  Rhino.59  Here,  the  axial  profile  is  rail  revolved  around 
the  base  cross  section  to  create  the  full  three-dimensional  surface.  This  surface  is 
exported  in  “.iges”  format  to  be  read  into  Pointwise  as  a  database.  This  database  is 
crucial  as  it  allows  for  grid  elements  to  be  resized  without  losing  important  geometric 
information.  A  surface  grid  is  created  using  this  database  as  a  reference,  and  it  is 
sized  and  broken  down  into  blocks  depending  on  the  needs  of  the  given  problem. 
Once  the  surface  grid  is  completed,  the  elliptical  solver  is  run  to  a  desired  smoothness 
using  whatever  control  functions  provide  the  best  solution  for  the  given  geometry. 
Finally,  a  volume  is  extruded  from  the  surface  hyperbolically,  with  user  prescribed 
boundary  conditions  and  spacing,  to  form  the  grid  to  be  exported  to  DPLR  (in 
Plot3D  format). 
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2.2.2  DPLR-specific  Grid  Generation  Concerns 


There  are  certain  known  issues  that  need  to  be  carefully  considered  when  craft¬ 
ing  a  blunt-body  topology  in  order  to  prevent  spurious  results  in  DPLR.  First,  a 
geometry  created  by  simply  sweeping  an  axial  profile  around  the  central  axis  of  a 
base  cross  section  will  produce  a  singularity  at  the  nose  of  the  resulting  vehicle.  This 
singularity  may  introduce  unwanted  errors  into  DPLR  (a  finite  volume  solver)  as 
cells  extruded  from  that  point  will  have  vanishing  volume.  In  order  to  remove  this 
singularity,  all  grids  are  patched  elliptically  in  the  nose  region  as  shown  in  Figure  2.4. 
Also,  in  order  to  prevent  poor  shock  capturing,  all  grids  use  hirer  spacing  in  highly 
curved  areas  (i.e.  shoulder  regions)  as  shown  in  figure  2.5.  Finally,  DPLR  requires 
a  sufficiently  small  body  normal  spacing  near  the  wall  in  order  for  the  hypersonic 
boundary  layer  to  be  captured.  To  this  end,  all  meshes  have  80  body  normal  points 
with  a  near  wall  spacing  of  1.0  x  10-6  nr. 
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Figure  2.4:  Nose  patching  as  shown  on  the  Apollo  heat-shield 


Figure  2.5:  Edge  spacing  for  the  Apollo  heat-shield 
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2.3  Shock  Layer  Radiation  Solver 

As  DPLR  lacks  a  method  for  internally  calculating  shock  layer  radiation,  an 
external  code  is  used  to  calculate  the  influence  of  the  radiating  flow  held  on  a  given 
a  blunt-body.  In  this  work,  the  Nonequilibrium  Air  Radiation  (NEQAIR)60  code 
is  used.  Radiation  is  only  calculated  in  an  uncoupled  sense,  in  that  results  from 
NEQAIR  are  not  fed  back  into  DPLR.  Essentially,  radiation  is  calculated  using 
only  the  fully  developed  solution  from  DPLR,  which  is  equivalent  to  applying  an 
optically  thin  assumption,  in  that  no  radiation  is  absorbed  by  the  shock  layer  itself. 
It  is  possible  to  loosely  couple  these  two  software  packages,  via  the  radiation  term  in 
the  conservation  of  energy  equations  as  in  the  work  by  Pace61  for  axisymmetric  test 
cases.  This  procedure  involves  iteratively  creating  new  radiation  solutions  based  on 
DPLR  solutions  updated  with  NEQAIR  data.  Full  three  dimensional  cases  using 
this  approach  can  be  incredibly  costly;  and,  since  the  work  of  Johnson  presents 
radiation  and  convection  in  an  uncoupled  sense,  it  is  reasonable,  for  the  sake  of  an 
apples-to-apples  comparison,  to  do  the  same  here. 

2.3.1  Radiation  Model 

NEQAIR  works  by  solving  the  radiative  transport  equation  (RTE),  given  by:38 

<^-=e-k'I  (2.19) 

as 

where  /  is  a  radiative  intensity,  e  and  k'  are  emission  and  absorption  coefficients 
respectively,  and  s  is  path  known  as  a  line-of- sight.  To  compute  these  coefficients, 
NEQAIR  uses  spontaneous  emission,  absorption,  and  stimulated  emissions  due  to 
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changes  in  energy  states  computed  along  a  linc-of-sight  for  every  chemical  species 
present.  Bound- free  and  free- free  radiation  is  also  considered.60  To  determine  the 
electronic  state  distribution,  the  quasi-steady  state  method  of  Park38  is  used,  which 
considers  electron  impact  excitation,  de-excitation,  and  recombination  in  forming 
a  model  for  the  population.  The  solution  of  the  radiative  transport  equation  is 
vastly  simplified,  within  NEQA1R,  by  applying  the  tangent-slab  approximation. 
This  model  essentially  makes  the  problem  one- dimensional  by  assuming  the  radiating 
shock  layer  to  be  an  infinitely  long  slab  of  radiating  gas  parallel  to  the  body  at  a 
specified  point.  As  such,  emission  and  absorption  can  be  neglected  in  the  body 
parallel  direction,  leaving  only  the  body  normal  (line  of  sight)  direction  upon  which 
to  integrate  the  intensities  in  the  RTE.  This  model  produces  heating  estimates  that 
are  5-15%  greater  that  what  would  be  predicted  by  models  that  include  surface 
curvature,62  with  the  advantage  of  significant  cost  savings. 

2.3.2  Work- flow 

A  NEQA1R  solution  is  generated  as  follows.  First,  species  concentrations  and 
temperatures  are  extracted  for  the  final  volume  mesh  of  a  converged  blunt-body  so¬ 
lution  using  the  built  in  post-processor  in  DPLR.  Then  lines-of-sight,  discretized  to 
a  set  number  points,  are  created  linearly  from  surface  nodes  to  the  outer  boundary 
of  the  volume  grid.  Thermodynamic  properties  are  interpolated  from  the  converged 
DPLR  solution  onto  corresponding  points  along  these  lines-of-sight.  Example  lines- 
of-sight  are  shown  in  Figure  2.6,  generated  for  the  Apollo  heat  shield.  Here,  2,667 
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lines  of  sight  are  needed  to  cover  the  part  of  the  vehicle’s  surface  needed  for  the 
simulation  (symmetry  dictates  that  only  half  of  the  heat  shield  is  needed).  Finally, 
NEQAIR  is  run  to  integrate  the  transport  equation,  along  each  line  of  sight,  to 
calculate  the  radiative  heat  flux  at  the  originating  surface  node  point  due  to  ra¬ 
diative  phenomena.  For  example,  the  topology  in  Figure  2.6  needs  2,667  separate 
invocations  of  NEQAIR  to  form  a  complete  solution. 


Figure  2.6:  Lines-of-sight  for  an  Apollo  heat- shied 
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2.4  Hardware 


The  above  computational  tools,  namely  DPLR  and  NEQAIR,  need  a  super¬ 
computing  environment  in  which  to  conduct  their  simulations.  Two  separate  plat¬ 
forms  are  used  primarily  in  this  work:  the  (1)  Skystreak  cluster  at  the  University  of 
Maryland  and  (2)  the  Columbia  cluster  at  NASA  Ames  Research  Center.  Skystreak 
,  upon  which  most  DPLR  simulations  are  run,  is  a  Gentoo  Linux  based  system  that 
consists  of  7  dual  processor  32-bit  AMD  Opteron  nodes  and  7  dual  processor  64- 
bit  Opteron  nodes.  This  setup  allows  for  a  maximum  14  computational  nodes  per 
problem;  and,  in  turn,  DPLR  solutions  require,  in  general,  840  to  4360  CPU  hours 
to  reach  final  convergence  (after  multiple  individual  runs  of  DPLR).  Skystreak  uses 
the  PathScale48  suite  for  Fortran  compiling  and  OpenMPI49  is  for  parallelization. 

The  Columbia  supercomputer,  used  for  NEQAIR  simulations  and  some  DPLR 
grid  resolution  cases,  is  capable  of  88.88  teraflops  per  seconds  using  13,312  total 
computational  cores  and  a  SUSE  Linux  operating  system.  It  is  made  of  17  Altix 
3700  (512  cores  each)  nodes  and  4  Altix  4700  nodes  (3584  total  cores)  nodes.  Intel 
Fortran50  and  OpenMPI  form  the  compilation  environment  on  Columbia.  Typical 
NEQAIR  runs  using  this  cluster  take  approximately  15  minutes  per  line  of  sight, 
while  DPLR  solutions  run  on  the  order  of  what  is  experienced  on  Skystreak  with  the 
bonus  that  more  processors  can  be  expended  on  a  given  problem. 
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Chapter  3 


Apollo  4  Benchmarking 

This  chapter  explores  benchmark  DPLR  solutions  using  the  Apollo  heat-shield 
at  Apollo  4  (AS-501)36,63  peak  heating  conditions  in  order  to  better  understand  po¬ 
tential  issues  that  may  be  encountered  with  DPLR  using  more  sophisticated  shapes. 
To  that  end,  the  baseline  Apollo  heat  shield  torus  is  altered  parametrically  to  fully 
comprehend  the  effects  of  certain  geometric  features,  particularly  edge  radius,  since 
shapes  generated  using  the  procedure  outlined  in  Chapter  1  do  not  posses  this  at¬ 
tribute. 

3.1  Baseline  Geometry  and  Design  Point 

The  baseline  Block  1  Apollo  command  module  geometry  is  shown  in  Figure 
3.1.  The  forebody  of  the  Apollo  command  module  during  re-entry,  or  the  portion 
of  the  vehicle  that  is  composed  of  the  heatshield,  is  a  23°  half-angle  spherical  seg¬ 
ment  blended  into  a  torus  with  radius  Rt  =  0.196  m  that  extends  for  133.9°.  The 
afterbody  is  a  33°  conical  frustum  with  a  cylindrical  cap  (for  the  purposes  of  this 
work,  the  outer  mold  line  of  the  Apollo  capsule  is  assumed  to  posses  a  spherical  cap 
with  radius  R  =  0.231  m).  For  basic  stability  comparisons,  the  center  of  gravity  is 
taken  at,  with  respect  to  the  nose  of  the  vehicle,  at  xcg  =  1.35  m,  ycg  =  0.00  m, 
zcg  =  —0.137  m,  consistent  with  what  was  used  for  wind-tunnel  tests.64 
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For  computational  purposes,  a  four-block  singularity  free  mesh  with  195,840 
grid  cells  (80  points  body  normal)  is  used  as  shown  in  Figure  3.2.  The  grid  consists 
of  the  Apollo  heat-shield  cut  off  at  its  widest  extent,  retaining  only  the  forebody  and 
neglecting  everything  in  the  afterbody.  Because  this  geometry,  and  all  subsequent 
geometries  considered,  is  symmetric  about  the  x-y  plane,  only  half  the  heat  shield 
needs  to  be  included  computationally,  thus  saving  on  computational  costs.  The 
baseline  DPLR  simulations  of  this  mesh  are  conducted  at  Apollo  4  peak  heating 
conditions,63  experienced  at  an  altitude  of  61  km,  M, ^  =  32.8,  and  a  =  —25°. 


50 


(a)  Surface 


(b)  Symmetry  Plane  (every  other  point) 


Figure  3.2:  Apollo  heat  shield  CFD  mesh 
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3.2  Baseline  DPLR  Results 


Figure  3.3  shows  Mach  contours  on  the  symmetry  plane  and  pressure  con¬ 
tours  on  the  surface  for  the  baseline  Apollo  4  case.  Lift  and  drag  coefficients  are 
calculated  as  0.465  and  1.220  respectively,  yielding  an  L/D  ratio  of  0.381  which  is 
very  close  to  the  ratio  of  0.375  predicted  by  flight  data65  at  this  particular  Mach 
number.  Using  the  center  of  mass  defined  above,  the  moment  coefficient,  Cm:Cg,  is 
0.016  which  matches  well  with  high  speed  wind  tunnel  data.64  Figure  3.4  shows 
the  total  wall  heat  flux  on  the  vehicle’s  surface  along  its  plane  of  symmetry.  Peak 
convective  heating  occurs  on  the  windward  (the  part  of  the  vehicle  pointing  into  the 
wind  when  at  angle  of  attack)  heat  shield  edge,  away  from  the  stagnation  point  with 
a  value  of  approximately  371  W/cm2,  1.68  times  the  stagnation  point  value.  This 
result  is  consistent  with  a  combination  of  the  observation  by  Lee  and  Goodrich36 
that  the  maximum  convective  heat  flux,  at  zero  angle  of  attack,  is  60%  larger  than 
at  the  stagnation  point  and  the  1.06  correction  suggested  by  Bertin8  to  account  for 
sonic  line  movement  when  the  vehicle  is  pitched.  Table  3.1  shows  a  comparison  of 
the  present  calculated  convective  heating  rates  with  past  work.  The  lower-order 
approach  underestimates  computational  solutions  by  30%  for  both  peak  and  stag¬ 
nation  point  convective  heat  flux.  These  under-predicted  values  by  the  lower-order 
method  compared  to  DPLR  solutions  can  be  attributed  to  the  failure  of  the  later  ap¬ 
proach  to  account  for  boundary  layer  blowing.  In  fact,  the  classic  Fay  and  Riddell15 
solution,  which  does  account  for  boundary  layer  physics,  as  calculated  by  Park66  is 
almost  identical  to  the  DPLR  result.  Still,  for  those  solutions  in  which  radiation  is 
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coupled  to  the  convective  heat  fluxes,  estimates  over-predict  DPLR  solutions  by  up 


to  65%. 
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(a)  Surface 


(b)  Symmetry  Plane 


Figure  3.3:  Pressure/convective  heating  contours  on  surface  and  Mach  contours  in 


symmetry  plane  for  Apollo  4  at  peak  heating  conditions 
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400 


Figure  3.4:  Symmetry  line  heating  profile  for  Apollo  4  at  peak  heating  conditions 


Table  3.1:  Apollo  4  peak  and  stagnation  point  convective  heating 


Author 

Qmax,conv{^^ /dll  ) 

%Diff 

qS,conv  (W/cm2) 

%Diff 

Present  work 

371 

221 

Johnson  et  al.2 

260 

-29.9 

154 

-30.3 

Pavlovsky  and  Leger67 

266“ 

-28.3 

- 

- 

Fay  and  Riddell15 

- 

- 

230 

4.1 

Park66 

- 

- 

363 

64.3 

Curry  and  Stephens68 

- 

- 

339 

53.4 

Bartlett  et  al.69 

- 

- 

289 

30.8 

Ried  et  al.21 

“ 

- 

227 

2.7 

“After  subtracting  qrad  from  Tauber  and  Sutton1' 
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At  high  entry  velocities,  shock  layer  radiation  becomes  significant  for  blunt 
bodies  clue  to  ionization  and  dissociation;69  and,  as  of  now,  DPLR  has  only  limited 
internal  methods  for  calculating  the  heating  impact  from  the  shock  layer.  Since 
the  lower-order  approach  does  not  couple  convective  and  radiative  heat  transfer, 
comparing  the  convective  heat  flux  is  an  informational  point  of  comparison  between 
the  two  aerothermodynamic  calculations,  yet,  if  accuracy  is  required,  the  present 
solutions  would  need  to  be  coupled  to  an  external  shock  layer  radiation  code.  The 
true  nature  of  the  coupling  between  radiative  and  convective  heat  transfer  (as  ev¬ 
idenced  by  the  spread  in  estimates  by  Park,  Bartlett  et  ah,  and  Reid  et  al.)  are 
not  well  understood  for  these  high  temperature  environments;  as  such,  it  is  reason¬ 
able  to  leave  the  two  heating  regimes  uncoupled  until,  at  least,  more  flight  data  is 
accrued  and  better  correlations  are  derived.  Essentially,  though  DPLR  convective 
heat  transfer  results  are  not  entirely  accurate  (in  that  shock  layer  radiation  is  ne¬ 
glected),  the  observation  that  peak  heating  is  not  at  the  stagnation  point,  and  is,  in 
fact,  significantly  higher  shows  that  computational  solution  yields  results  that  are, 
in  the  least,  qualitatively  significant. 

3.3  Grid  Resolution 

Solutions  for  the  Apollo  heat  shield  at  Apollo  4  (AS-501)  peak  heating  condi¬ 
tions,  experienced  at  an  altitude  of  61  km,  M, x  =  32.8  at  a  =  —25°,  are  compiled 
on  volume  meshes  with  40,  60,  80,  and  160  points  in  the  body  normal  direction. 
The  80  point  solution  is  the  baseline  case  used  in  benchmark  DPLR  runs.  Table  3.2 
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shows  point  aerothermodynamic  metrics  for  all  these  cases,  and  percent  differences 
are  referenced  to  the  highest  resolution  case  (160  points).  The  lowest  resolution  case 
shows  the  greatest  disparity  with  respect  to  the  finer  case,  especially,  with  regards 
to  peak  convective  heat  flux.  Here  the  difference  is  greater  than  10%  whereas  adding 
just  20  more  body  normal  points  drops  the  disparity  below  1%.  Both  the  80  point 
case  and  60  point  case  display  errors  less  than  1%  for  all  aerothermodynamic  metrics 
compared  with  the  finest  resolution  case;  however,  not  all  heat  shields  have  a  simple 
axisymmetric  shape  like  Apollo,  nor  do  these  point  metrics  accurately  portray  what 
is  happening  on  the  entire  surface. 

Figure  3.5  shows  the  convective  heat  flux  on  the  symmetry  plane  for  all  four 
grid  resolution  cases.  Clearly,  the  lowest  resolution  case  shows  the  most  erroneous 
predictions,  especially  on  the  leeward  side  of  the  vehicle  (the  portion  of  the  heat 
shield  that  points  away  from  the  wind  when  pitched).  One  interesting  note  is  that  all 
solutions,  even  the  sparsest  mesh  case,  converge  at  the  stagnation  point  convective 
heat  flux.  The  60  and  80  point  cases  are  nearly  identical,  but  with  the  hirer  solution 
portraying  slightly  better  results  on  the  leeward  side  of  the  vehicle.  In  this  instance, 
results  reported  from  an  Apollo  mesh  with  60  points  in  the  body  normal  would  be 
nearly  identical  to  the  baseline;  however,  since  not  all  heat  shields  studied  here  are 
simple  shapes,  a  slightly  larger  resolution  is  a  safer  choice.  In  that  regard,  choosing 
80  points  for  the  body  normal  direction  for  all  heat  shield  meshes  is  a  prudent  one. 
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Table  3.2:  Apollo  grid  resolution  aerothrmodynamics 


Parameter 

40 

Pts  (body 
60 

normal) 

80* 

160 

CL 

0.457 

0.463 

0.463 

0.467 

A(%) 

-2.141 

0.463 

0.463 

- 

C D 

1.20 

1.22 

1.22 

1.23 

A(%) 

-2.44 

-0.81 

-0.81 

L/D 

0.3808 

0.3795 

0.3795 

0.3797 

A(%) 

0.31 

-0.04 

-0.04 

- 

r 

w m,cg 

0.0159 

0.0163 

0.0162 

0.0161 

A(%) 

-1.24 

1.24 

0.62 

- 

Qconv,max  ('W'/cm  ) 

416.5 

373.4 

373.8 

370.5 

A(%) 

12.42 

0.78 

0.89 

“ 

^Baseline  dimension  for  all  DPLR  solutions 


y  (m) 


Figure  3.5:  Symmetry  plane  convective  heat  flux  for  Apollo  heat  shield  grid  resolu¬ 


tion  study 
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3.4  Baseline  Radiation  Results 


NEQAIR  surface  radiative  heating  results  generated  from  the  baseline  DPLR 
solution  are  shown  in  Figure  3.6.  Peak  radiative  heating  occurs  slightly  leeward  of 
the  stagnation  point  at  a  value  of  217  W/cm2.  Stagnation  point  radiative  heating 
is  only  slightly  lower  at  204  W/cm2.  Table  3.3  shows  a  comparison  of  the  present 
calculated  stagnation  point  radiative  heating  rates  with  past  work.  The  lower-order 
method  over  estimates  NEQAIR  stagnation  point  radiative  heat  flux  by  16%;  how¬ 
ever,  this  overestimation  combined  with  the  under-predicted  convective  heat  flux 
yields  a  total  stagnation  point  heat  flux  (391  W/cm2)  that  differs  from  the  com¬ 
bined  DPLR/NEQAIR  solution  (425  W/cm2)  by  only  8%.  Total  stagnation  point 
heat  flux  agrees  mostly  well  (within  25%)  with  past  results,  while  stagnation  point 
radiative  heat  flux  agrees  to  within  20%  for  most  cases.  The  values  that  are  vastly 
greater  than  the  NEQAIR  solution  use  models  that  do  not  account  for  energy  dis¬ 
sipation  in  the  boundary  layer,  resulting  in  expectedly  greater  estimates.  It  should 
be  noted  that  for,  convective  heating,  past  work  resulted  in  higher  values,  while 
that  trend  is  the  opposite  for  radiative  heating.  Combining  these  two  phenomena 
explains  the  relative  accuracy  of  the  combined  DPLR/NEQAIR  approach  to  predict 
total  stagnation  point  heat  flux.  This  implies  that  the  effect  of  coupling  would  be 
to  lower  radiative  heat  flux,  while  increasing  convective  heat  flux,  creating  only  a 
negligible  difference  in  the  sum.  Figure  3.7  shows  the  convective,  radiative,  and 
total  heat  flux  along  the  symmetry  plane.  The  maximum  overall  heat  flux  occurs 
at  the  location  of  peak  convective  heating  (windward  of  the  stagnation  point)  at  a 


value  507  W/cm2.  This  prediction  agrees  very  well  with  the  480  W/cm2  estimation 


by  Pavlovsky  and  Leger.67 
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Figure  3.6:  Surface  radiative  heat  flux  contours  for  Apollo  heat  shield  at  Apollo  4 


peak  heating  conditions 
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Tabic  3.3:  Apollo  4  stagnation  point  radiative  and  total  heating 


Author 

qS,rad(  W/cm2) 

%Diff 

qs,tot(  W/cm2) 

%Diff 

Present  work 

204 

425 

Johnson  et  al.2 

237 

16.18 

391 

-8.00 

Flight  data0. 

167 

-18.1373 

- 

- 

Curry  and  Stephens68 

176 

-13.73 

515 

21.18 

Bartlett  et  al.69 

193 

-5.39 

482 

13.41 

Ried  et  al.21 

300 

47.06 

527 

24.00 

Park  200 170 

507 

148.53 

731 

72.00 

Park  200466 

168 

-17.65 

531 

24.94 

LORAN71 

320 

56.86 

- 

- 

Balakrishnan  et  ah' 2 

184 

-9.80 

- 

“  Calculated  using  measured  peak  intensity  by  Park76 


Figure  3.7:  Convective,  radiative,  and  total  heat  flux  on  symmetry  plane  for  Apollo 
heat  shield  at  Apollo  4  peak  heating  conditions 
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3.5  Torus  Radius 


As  previously  stated,  the  forebody  of  the  Apollo  command  module  during  re¬ 
entry  is  made  up  of  a  23°  half-angle  spherical  segment  blended  into  a  torus  with 
radius  Rt  =  0.196  m  that  extends  for  133.9°  until  the  conical  frustum  afterbody 
begins.  This  torus,  only  to  the  capsules  widest  extent,  is  included  in  the  baseline 
benchmark  results  for  Apollo  4;  however,  the  geometries  generated  by  the  optimiza¬ 
tion  process  do  not  include  these  regions  of  curvature.  What  role  this  torus  plays  on 
a  blunt-body’s  flow  held  is  detailed  in  the  following  subsections  through  parametric 
studies  of  the  baseline  Apollo  geometry  with  various  different  torus  designs. 

3.5.1  Torus  Extent 

To  further  understand  the  effects  the  extent  of  this  torus  has  on  the  aerother- 
modynamics  of  a  blunt-body  capsule,  the  baseline  Apollo  4  case  (referred  to  here 
as  a  half  torus)  is  compared  to  solutions  on  meshes  that  include  the  entire  torus 
(expanded  until  the  beginning  of  the  conical  frustum  afterbody)  and  that  include 
no  torus  at  all.  Figure  3.8  shows  a  comparison  of  the  wall  heat  flux  on  the  plane 
of  symmetry  for  all  three  cases.  Only  minor  differences  can  be  seen  when  any  part 
of  the  torus  is  considered;  however,  the  results  obtained  without  a  torus  show  a 
singularity  at  the  edges  of  the  heat  shield  and  a  local  maximum  in  a  very  different 
location  than  the  other  cases.  Since  local  convective  heating  is  proportional  to  the 
reciprocal  of  the  square  root  of  the  local  radius  of  curvature,  this  asymptotically 
high  heating  at  the  edge  of  the  heat  shield  is  not  surprising;  however,  another  ex- 
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planation  for  these  spurious  results  can  be  seen  by  examining  the  sonic  line  of  the 
three  cases  (Figure  3.9). 

For  cases  that  include  the  torus,  regions  of  subsonic  flow  spill  over  to  the  torus 
on  the  windward  side  of  the  heat  shield  before  expanding  back  to  supersonic.  When 
the  torus  is  absent,  flow  at  the  exit  of  the  grid  near  the  windward  edge  is  subsonic, 
violating  the  DPLR  boundary  condition  of  a  supersonic  exit.  This  numerical  limita¬ 
tion  can  be  solved,  without  introducing  some  curvature  at  the  edges,  by  adding  an 
afterbody  to  the  heat  shield.  Still,  by  adding  a  simple,  conical  afterbody  (thereby 
retaining  the  infinitesimally  small  radius  of  curvature  at  the  heat  shield  edge)  with¬ 
out  including  the  wake,  there  is  still  the  possibility  of  subsonic  flow  at  the  leeward 
boundaries  of  the  heat  shield.  Essentially,  this  means  that  the  afterbody  must  be 
a  fully  closed  body  and  the  topology  must  be  extended  to  include  the  wake,  which 
will  greatly  increase  the  computational  cost  incurred.  To  avoid  these  sky-rocketing 
costs,  some  manner  of  curvature  at  the  edge  of  an  optimized  heat  shield  must  be 
added  while  taking  care  to  not  drastically  change  the  original  geometry. 
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y  (m) 


Figure  3.8:  Symmetry  convective  heat  flux  comparison  of  Apollo  4  peak  heating 


case  for  three  different  torus  extents 


Figure  3.9:  Sonic  line  comparison  of  Apollo  4  peak  heating  case  for  three  different 


torus  extents 
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3.5.2  Toms  Size 


Though  some  finite  curvature  is  required  to  perform  CFD  with  DPLR  on 
a  blunt-body  heat  shield,  it  was  unclear  what  impact  the  exact  amount  of  edge 
curvature  or  bluntness  will  have  on  the  aerothermodyamics  of  a  heat  shield.  As 
such,  Apollo  forebodies  with  various  multiples  of  the  baseline  torus  radius  (up  to 
5 xRt)  are  compared  to  the  benchmark  Apollo  results  calculated  using  Apollo  4 
peak  heating  altitude  and  velocity  at  the  dominant  entry  angle  of  attack,  a  =  —25°, 
and  at  a  =  —15°.  The  lower  magnitude  angle  of  attack  case  is  included  to  explore 
edge  radius  effects  on  a  blunt-body  for  which  the  stagnation  point  is  not  near  the 
vehicle’s  edge.  Figure  3.10  shows  the  maximum  heat  flux  on  the  symmetry  plane  for 
all  topologies  at  both  angles  of  attack.  Both  curves  show  that  heat  flux  decreases 
in  a  power  law  sense  with  increasing  torus  radius,  which  should  be  the  case  since 
convective  heat  transfer  is  function  of  the  inverse  of  the  square  of  the  local  radius  of 
curvature.  However,  when  the  stagnation  point  is  further  from  the  vehicles  windward 
edge  (i.e.  for  a  =  —15°)  the  heat  flux  decreases  at  a  slightly  slower  rate  than  it  does 
at  a  higher  angle  of  attack.  Since  the  regions  of  higher  temperature  and  pressure 
occur  farther  away  from  the  torus  at  lower  angles  of  attack,  it  is  reasonable  to 
assume  that  changing  the  torus  radius  will  have  a  lesser  impact  on  the  heat  flux 
in  this  case.  For  all  torus  sizes,  the  a  =  —15°  case  has  lower  peak  convective  heat 
flux.  This  phenomena  can  be  attributed  to  the  lower  velocity  gradients  experienced 
at  the  edge  of  the  heat  shield  at  that  particular  angle  of  attack. 

Figure  3.11  shows  the  resulting  lift  to  drag  ratio  for  all  cases  of  torus  radius  for 
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both  considered  angles  of  attack,  -15°and  -25°.  The  dashed  lines  show  the  predicted 
Newtonian  L/D  for  each  respective  angle  of  attack.  When  torus  radius  is  small,  the 
DPLR  calculated  L/D  approaches  this  predicted  value.  Since  the  Newtonian  values 
are  calculated  with  no  radius  of  curvature  at  the  edge,  this  trend  is  understandable. 
For  both  angles  of  attack,  the  L/D  decreases  linearly;  however  the  lift  to  drag  ratio 
for  the  higher  angle  of  attack  set  of  solutions  decreases  at  a  faster  rate.  This  occurs 
due  to  peak  pressures  shifting  toward  the  center  of  the  heat  shield  rather  than  being 
closer  to  the  highly  curved  edges.  Figure  3.12  is  a  plot  of  moment  coefficient  versus 
torus  size  for  both  considered  angles  of  attack.  The  center  of  gravity  may  change 
with  increasing  torus  radius  if  more  heat  shield  material  is  added  unevenly  to  the 
windward  edge  to  counter  higher  heating;  but,  for  all  cases  studied  here,  the  center 
of  gravity  is  considered  fixed.  Cmfig  increases  linearly  with  torus  radius  for  both 
angles  of  attack  due  to  the  increased  moment  arm  induced  by  the  larger  tori.  Also, 
note,  that  for  a  larger  edge  radius  (~  5 xRt),  ol  =  —15°  is  nearly  the  trim  angle  of 
attack. 
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Figure  3.10:  Peak  convective  heat  fluxes  for  Apollo  heat  shield  at  Apollo  4  peak 
heating  conditions  for  two  angles  of  attack 


Figure  3.11:  Lift  to  drag  ratios  versus  torus  radius  for  Apollo  heat  shield  at  two 
angles  of  attack 
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Figure  3.12:  Moment  coefficient  for  Apollo  heat  shield  at  different  angles  of  attack 

3.5.3  Further  Considerations 

Edge  curvature  plays  a  significant  role  in  determining  the  resulting  aerother- 
modynamics  of  a  vehicle.  Mission  requirements  and  hardware  concerns  (i.e.  launch 
vehicle  mating)  will  often  determine  what  the  afterbody  of  a  vehicle  will  look  like 
and  how  it  will  attach  to  the  heat  shield.  Since  these  concerns  are  beyond  the  scope 
of  this  work,  a  fixed  torus  matching  that  of  the  Apollo  capsule  (Rt  =  0.196  m)  is 
blended  to  all  blunt-body  optimized  designs  for  further  study.  However,  for  elliptical 
bases,  scaling  effects  need  to  be  considered  when  deciding  upon  which  axis  to  apply 
the  torus.  Figure  3.13  shows  an  elliptical  heat  shield  (to  be  discussed  in  Chapter  4) 
with  the  torus  applied  on  the  semi-major  axis  on  the  left  and  on  the  semi-minor  axis 
on  the  right.  Although  applying  the  curvature  along  the  shorter  side  first  generates 
a  larger  surface  area,  the  peak  heat  flux  is  significantly  less  than  if  the  torus  were 
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first  applied  on  the  longer  side.  These  lower  heating  rates  would,  in  all  likelihood, 


require  cheaper  materials  to  dissipate.  As  such,  all  elliptical  and  blunt  designs  to 
be  studied  here  will  have  the  torus  added  to  the  axial  profile  along  the  semi-minor 
axis  before  being  swept  around  the  base  cross-section. 


Original  shape 


q.  (W/cm2):  20  80  140  200  260  320  380 


Figure  3.13:  Convective  heat  flux  on  an  elliptical  heat  shield  for  two  different  meth¬ 
ods  of  torus  generation 

3.6  Computational  Cost  Summary 

Table  3.4  shows  a  summary  of  grid  sizes,  iteration  counts,  and  computational 
time  for  all  Apollo  derived  cases  considered  here.  There  is  no  consistency  in  the 
time  it  takes  to  arrive  at  a  final  solution,  since  convergence  occurs  only  at  the  user’s 
satisfaction.  That  is,  the  solution  is  allowed  to  mature  until  the  user  deems  the 
solution  to  have  stabilized  (usually  by  the  time  the  L2  norm  of  the  residuals  of  a 
conserved  variables  is  less  than  10~8).  Still,  no  case  took  less  than  390  CPU  hours 
and  20,000  iterations  for  convergence.  In  fact,  most  cases  needed  more  time  to  reach 
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sufficient  convergence.  NEQAIR  solutions  were  run  using  2,667  lines  of  sight  for  the 


Apollo  geometry.  At  15  minutes  per  line  of  sight,  the  single  uncoupled  radiation 
solution  needed  666.75  CPU  hours  to  complete. 


Table  3.4:  Summary  of  computational  cost  for  Apollo  edge  radius  cases 


R/RT 

a 

Grid  Size  (#  cells) 

Iterations 

CPU  Time  (hrs) 

0.100 

-25° 

197120 

64400 

1831.2 

0.200 

-25° 

197120 

77900 

1594.46 

0.200 

-15° 

197120 

52500 

987.84 

0.250 

-25° 

191720 

70000 

1470.0 

0.333 

-25° 

191720 

46600 

1470.0 

0.500 

-25° 

203040 

66000 

1186.08 

0.750 

-25° 

203040 

85300 

2286.62 

1.000 

-25° 

195840 

20400 

392.84 

1.000“ 

-25° 

98560 

108000 

1099.99 

1.000“ 

-25° 

147840 

37000 

604.33 

1.000“ 

-25° 

394240 

55500 

2123.33 

1.000 

-15° 

195840 

42900 

1524.46 

1.500 

-25° 

205560 

36800 

735.0 

1.500 

-15° 

205560 

34000 

692.16 

2.000 

-25° 

205560 

41700 

976.08 

2.000 

-15° 

205560 

40700 

1240.54 

2.500 

-25° 

205560 

38100 

1034.46 

2.500 

-15° 

205560 

33600 

738.92 

5.000 

-25° 

205560 

52300 

1664.46 

5.000 

-15° 

205560 

31400 

1065.54 

“Grid  resolution  cases 
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Chapter  4 


Slender  Bodies:  A  High  L/D  Case 

The  n2  parameter,  from  Equation  1.1,  controls  the  sharpness  of  the  base  cross 
section  for  blunt-body  heat  shields  of  Johsnon,  et  al.4'6  Optimized  geometries  have 
utilized  a  lower  bound  of  1.3  for  this  parameter  along  with  a  more  slender  profile 
to  generate  high  L/D  solutions;  however,  it  is  unclear  as  to  how  this  sharpness  will 
affect  the  off-stagnation  point  performance  of  these  generated  shapes.  As  such,  a 
representative  geometry,  optimized  for  high  L/D  at  Apollo  4  peak  heating  condi¬ 
tions,  is  explored  in  this  chapter. 

4.1  Baseline  Geometry  and  Results 

To  study  the  full  flow  field  of  the  high  L/D  shapes  classified  by  Johnson,  et 
al.2,6  ,  DPLR  solutions  are  obtained  for  a  spherical  segment  heat  shield  optimized 
for  maximum  L/D  at  Apollo  4  peak  heating  conditions  ( ht  =  61  km,  =  32.8). 
The  modified  Newtonian  approach  predicted  a  lift  to  drag  ratio  of  1.24  and  qs,Conv  of 
240  W/cm2  for  an  89°  spherical  segment  with  n2  =  1.3,  m  =  4,  and  e  =  —0.968  at 
a  =  18°.  This  shape  (see  Figure  4.1)  has  more  of  a  “nosecone”  like  geometry  in  that 
it  is  more  slender  and  eccentric  than  the  Apollo  heat  shield.  As  such,  edge  radius 
effects  would  be  expected  to  have  less  significant  an  impact  on  overall  performance, 
since  the  velocity  gradients  at  the  boundaries  are  less  steep  than  those  for  a  more 
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classical  spherical  segment  blunt  design,  for  this  study  a  four  block  structured  mesh 


with  236,440  grid  cells  is  used  as  shown  in  Figure  4.2. 
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(a)  Surface 


(b)  Symmetry  Plane  (every  other  point) 


Figure  4.2:  mesh  for  89°  spherical  segment  optimized  for  maximum  L/ D 
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Figures  4.3  and  4.4  show  surface  heat  flux  and  pressure  contours  as  well  as 
symmetry  plane  Mach  and  pressure  contours  respectively.  Mach  contours  are  very 
close  to  the  vehicle’s  surface,  resulting  in  a  thinner  shock  layer  for  radiation  than 
exists  with  a  more  blunt  geometry  (i.e.  Apollo).  DPLR  solutions  show  Cl  =  0.795 
and  Cr,  =  0.711,  resulting  in  L/D  =  1.118.  The  lower-order  prediction  for  L/D 
is  11%  higher  than  the  CFD  result.  This  difference  is  due,  in  most  part,  to  the 
simpler  method  ignoring  surface  pressures  in  the  shadow  region  (V  ■  n  >  0)  that 
may  contribute  to  lift  and  drag.  Still,  this  offset  is  not  extreme;  and,  the  predictions 
by  the  lower-order  approach  would  still  be  useful  in  design  studies. 

Convective  heating  at  the  stagnation  point  is  calculated  as  430  W/cm2.  The 
low-order  approach  under-predicts  this  value  by  44%,  a  larger  discrepancy  than  was 
observed  for  Apollo  4.  This  difference  implies  that  empirical  correlations  used  for 
heating  rates  may  not  be  ideally  suited  for  elliptical  base  cross  sections  such  as  the 
one  seen  here.  Peak  heating  occurs  along  the  leading  edge  of  the  vehicle  far  away 
from  the  stagnation  point  at  a  value  of  approximately  980  W/cm2.  This  value  is 
2.28  times  higher  than  the  DPLR  calculated  stagnation  point  heat  transfer  and  4.08 
times  greater  than  the  low-order  stagnation  point  prediction.  Essentially,  the  edge 
of  the  parallelogram  base  cross  section  (controlled  by  the  77,2  parameter)  creates  a 
sharp  leading  edge  away  from  the  nose  near  the  edges  of  the  heat  shield,  that  gen¬ 
erates  what  is  almost  an  attached  shock-wave  at  the  point  where  highest  convective 
heating  is  shown  to  occur.  The  heat  flux  is  expectedly  high  in  that  area  because 
there  is  little  gas,  in  the  shock  layer,  with  which  to  dissipate  the  high  temperatures 
created  by  the  shock-wave.  Because  the  lower-order  approach  only  looks  at  the  stag- 
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nation  point  for  heating  rates,  it  omits  an  area  of  adverse  heating  that  may  make 
the  design  infeasible.  For  future  optimizations,  a  situation  such  as  this  one  can  be 
avoided  by  altering  the  optimization  constraints  in  such  a  way  that  would  prevent 
near  attached  shock- waves  from  forming  on  a  generated  heat  shield. The  simplest 
way  to  do  this  would  be  to  raise  the  lower  bound  of  the  sharpness  parameter,  n2  to 
1.4  or  1.5. 


segment  optimized  for  maximum  L/D 
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Figure  4.4:  Symmetry  plane  contours  for  89°  spherical  segment  optimized  for  max¬ 


imum  L/D 


75 


4.2  Surface  Grid  Resolution 


Furthermore,  the  impact  of  surface  grid  resolution  on  the  resulting  aerother- 
modynamics  of  the  above  slender  heat  shield  is  examined  by  comparing  results  from 
the  representative  high  L/ D  shape  (an  89°  spherical  segment  with  n2  =  1.3,  m  =  4, 
and  e  =  —0.968  at  a  =  18°  at  Apollo  4  peak  heating  freestream  conditions)  using 
a  baseline  surface  grid  (236,440  total  cells)  with  one  in  which  the  number  of  points 
in  both  surface  directions  (i  and  j)  are  doubled  (1,205,600  total  cells).  Table  4.1 
shows  the  resulting  aerothermodynamics  of  these  two  cases.  The  lower  resolution 
case  shows  good  agreement  (within  1.5%)  with  the  hirer  resolution  solution  for  aero¬ 
dynamic  coefficients  of  lift  and  drag  but  shows  a  greater  than  15%  difference  with 
respect  to  peak  convective  heating. 

A  closer  examination  of  the  surface  convective  heat  flux  of  these  two  meshes 
(see  Figure  4.5)  shows  a  clearer  picture  of  this  disparity.  In  both  cases,  peak  heating 
occurs  at  nearly  the  same  position,  on  the  leading  edge  away  from  the  stagnation 
point,  but  the  higher  resolution  solution  shows  a  larger  area  of  high  heating  around 
that  point  at  values  much  higher  than  seen  in  the  low  resolution  case.  Likely,  the 
greater  number  of  points  generates  a  shape  that  is  sharper  than  geometrically  possi¬ 
ble  with  fewer  points.  In  reality,  a  heat  shield  with  this  design  would  have  a  leading 
edge  that  would  ablate,  or  burn  up,  as  it  experienced  these  high  heat  loads.  Es¬ 
sentially,  it  would  have  an  initial  shape  that  is  more  like  the  fine  mesh  that  would 
eventually  become  more  like  the  sparser  mesh  over  time.  The  lower  point  case  is 
probably  more  artificially  blunt  than  it  is  meant  to  be,  if  only  the  geometric  parame- 
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ters  were  considered,  which  causes  heating  rates  to  be  reported  lower.  Still,  nothing 
in  these  results  conflicts  with  the  above  conclusion  that  an  ri2  value  of  1.3  allows 
for  the  possibility  of  extreme  off-stagnation  point  heating  rates  that  may  reach  or 
exceed  material  limits;  in  fact,  they  only  underscore  it. 


Table  4.1:  Aerothermo dynamics  for  slender  heat  shield  surface  grid  resolution  study 


Cells  (Surface) 

CL 

cD 

L/D 

Qconv,max  ("W / CHI  ) 

236,440 

0.795 

0.711 

1.118 

978.75 

1,205,600 

0.800 

0.705 

1.135 

1173.2 

A(%) 

-0.625 

0.851 

-1.464 

-16.574 

q™w  (W/cm2):  100  300  500  700  900  1100 


Figure  4.5:  Surface  convective  heat  flux  for  slender  heat  shield  surface  grid  resolution 
study 
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4.3  Radiation 


NEQAIR  surface  radiative  heating  results  generated  from  the  high  L/D  test 
case  DPLR  solution  are  shown  in  Figure  4.6.  Peak  radiative  heating  is  106  W/cm2 
and  occurs  at  the  stagnation  point.  The  lower-order  approach  predicts  stagnation 
point  radiative  heating  as  190  W/cm2,  which  is  79%  higher  than  the  NEQAIR 
result.  This  disparity  is  much  greater  than  the  16%  difference  seen  in  the  Apollo 
benchmarking  case.  As  such,  it  is  entirely  possible  that  the  assumptions  used  in 
applying  the  Tauber  and  Sutton  model  for  radiative  heating  in  the  lower-order 
method  may  be  incorrect.  Likely,  the  relations  used  to  calculation  shock  stand  off 
distance,  the  most  important  factor  in  determining  the  radiative  heating,  do  not 
account  for  elliptical  geometries  such  as  those  with  oblate  stretching  in  the  base 
cross-section  seen  here.  Since  these  models  are  empirical  in  nature,  this  implies  the 
necessity  for  further  wind  tunnel  and  flight  tests  to  provide  the  data  points  needed 
in  creating  an  approximation  that  has  even  greater  physical  basis. 

Using  the  value  of  radiative  heating  calculated  by  NEQAIR,  the  total  stagna¬ 
tion  point  heating  is  536  W/cm2.  This  result  compares  favorably  to  the  lower-order 
prediction  of  430  W/cm2  (a  -20%  difference).  In  reality,  the  convective  and  radia¬ 
tive  heating  are  drastically  under-predicted  and  over-predicted  by  the  lower-order 
method  respectively  when  compared  to  the  present  computational  approach.  So, 
in  essence,  the  errors  introduced  to  the  convective  heat  flux  model  by  the  elliptical 
cross  section  are  canceled  out  by  the  errors  in  predicting  true  shock  stand  off  dis¬ 
tance  by  the  radiative  heating  model.  While  certainly  not  intended,  this  inaccuracy 
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in  the  underlying  models  generates  overall  results  that  are  indeed  accurate,  and  can 
be  used  effectively  for  first-pass  design  studies. 

Not  surprisingly,  NEQAIR  predicts  little  to  no  radiative  heating  on  the  lead¬ 
ing  edge  of  the  vehicle  at  the  location  where  maximum  convective  heat  flux  occurs. 
This  result  confirms  the  earlier  assertion  that  the  shock-wave  must  be  very  close  to 
the  body  at  that  point.  Essentially,  because  the  shock  layer  is  so  small,  there  is  very 
little  high  temperature  gas  necessary  to  radiate  heat  back  to  the  vehicle’s  surface. 
Still,  the  peak  convective  flux  is  almost  twice  the  total  heating  felt  at  the  stagnation 
point  and  can  not  be  ignored,  as  it  nears  design  limits  for  the  Orion  CEV  capsule. 


qrad  (W/cm2):  10  40  70  100 


Figure  4.6:  Surface  radiative  heat  flux  and  pressure  contours  for  89°  spherical  seg¬ 
ment  optimized  for  maximum  L/D 
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4.4  Changing  Edge  Sharpness 


The  effective  edge  sharpness  of  the  high  L/D  example  is  further  studied  by 
comparing  the  baseline  (77-2  =  1.3)  case  to  other  designs  with  the  all  same  parameters 
except  a  different  value  of  n2.  Figure  4.7  shows  the  peak  heat  fluxes  for  four  different 
values  of  n2,  and  it  can  be  seen  that  heating  decreases  approximately  with  the  inverse 
cube  of  the  sharpness  parameter.  Essentially,  by  increasing  the  bluntness  of  the 
edge  on  the  base  cross  section,  large  reductions  in  peak  heat  flux,  when  compared 
to  the  baseline  (n2  =  1.3),  are  obtained.  Figure  4.8  shows  the  lift  to  drag  ratios 
for  the  four  cases  of  n2  studied  here.  In  all  cases,  the  modified  Newtonian  solution 
over-predicts  the  calculated  L/D  by  up  to  10%.  This  discrepancy  is  caused  by  the 
lower-order  method’s  consistent  inability  to  capture  reductions  in  lift  do  to  pressures 
experienced  on  the  body  in  the  vehicle’s  shadow  region,  which  are  neglected  by  the 
modified  Newtonian  approach.  This  phenomenon  is  apparent  in  that  the  modified 
Newtonian  approach  predicts  lift  coefficients  that  are  2-10%  more  than  the  high 
fidelity  model;  however,  the  accuracy  of  the  Modified  Newtonian  approach  does 
improve  as  the  base  cross  section  becomes  more  elliptical  (n2  «  2.0).  When  the 
base  is  more  like  a  parallelogram  (n2  ~  1. 1-1.3),  the  revolved  surface  will  have  a 
sharp  leading  edge  blending  into  a  blunt  nose.  A  Newtonian  solution  is  not  as  well 
suited  for  these  sharp  leading  edges;9  as  such,  the  blunter  edged  solutions  (n2  ~  2.0) 
should  be  more  accurate. 
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Figure  4.7:  Maximum  convective  heat  fluxes  for  89°  spherical  segment  with  varying 
712  parameter 
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Figure  4.8:  Lift  to  drag  ratios  for  89°  spherical  segment  with  varying  77,2  parameter 
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4.5  Computational  Cost  Summary 


Table  4.2  shows  a  summary  of  grid  sizes,  iteration  counts,  and  computational 
time  for  all  cases  considered  in  this  section.  In  order  to  properly  mesh  the  paral¬ 
lelogram  cross-section  geometry,  more  points  are  needed  than  the  simpler  Apollo 
heat  shield.  As  such,  more  computational  effort  is  needed  in  general  for  these  cases. 
No  case  took  less  than  900  CPU  hours  and  34,000  iterations  to  converge.  The  grid 
resolution  case  took  over  204  days  of  CPU  time  to  complete.  This  was  facilitated 
by  the  much  larger  Columbia  supercomputing  cluster,  as  the  use  of  48  processors 
in  parallel  lessened  the  physical  duration  to  eight  and  half  days.  All  other  cases 
were  performed  on  the  Skystreak  cluster,  using  only  14  processors  in  parallel.  The 
NEQAIR  radiation  case  require  3,864  lines  of  sight  and  966  CPU  hours  to  complete. 


Table  4.2:  Summary  of  costs  for  89°spherical  segment  optimized  for  max  L/D 


ri2 

Grid  Size  (#  cells) 

Iterations 

CPU  Time  (hrs) 

1.1 

236440 

56000 

1625.56 

1.3 

236440 

65300 

1738.38 

1.3“ 

1205600 

86200 

9800.00 

1.5 

236440 

42000 

1135.54 

1.7 

236440 

52800 

1493.38 

2.0 

236440 

34500 

937.16 

“Grid  resolution  case  done  on  Columbia  cluster 
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Chapter  5 


Vehicle/Trajectory  Optimized  Geometries 

5.1  Lunar  Return  Optimized  Designs 

This  section  presents  computational  solutions  for  heat  shields  that  were  gen¬ 
erated  for  earth  return  following  a  mission  to  the  moon  using  the  trajectory/vehicle 
coupled  optimization  scheme  discussed  in  Chapter  1.  Due  to  time  and  computa¬ 
tional  constraints,  no  NEQA1R  radiation  simulations  are  undertaken  for  these  cases. 
As  such,  only  convective  heating  and  aerodynamic  calculations  are  presented. 

5.1.1  General  Summary 

CFD  solutions  are  obtained  for  heat  shields  generated  using  the  coupled  opti¬ 
mization  technique  for  lunar  return  entry  velocities,  Ve  =  11  krn/s,  and  entry  flight 
path  angles  of  7 e=  —6.0°  at  the  location  on  the  trajectory  upon  which  the  peak 
instantaneous  heat  flux  is  predicted  to  occur.  Table  5.1  shows  a  summary  of  the 
geometry,  design  point,  the  aerothermodynamics  calculated  from  DPLR  solutions 
and  the  predicted  aerothermodynamics  using  the  lower-order  approach  for  all  cases 
studied  in  this  section.  The  cases  maintain  their  descriptors  from  Table  12.1  of  Ref¬ 
erence  [6].  Percent  differences,  in  reference  to  the  DPLR  solutions,  are  presented  for 
the  lower-order  results  in  parenthesis  where  applicable.  All  shapes  considered  here 
have  a  spherical  segment  axial  profile,  as  all  other  choices  for  axial  shape  mimicked 
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a  spherical  segment  as  they  approached  an  optimum.  Basically,  optimized  power- 
law  and  spherically  blunted  cone  geometries  were  disguised  spherical  segments  for 
lunar  return  trajectories.  Optimizations  were  performed  at  L/ D  =  0.3,  0.5,  and  1.0 
for  two  objective  functions  sets:  (1)  maximizing  downrange  while  minimizing  heat 
load  and  (2)  maximizing  cross  range  while  minimizing  heat  load.  Trajectory  entry 
corridors  widths  of  up  to  1.37°  were  used  to  ensure  mission  feasibility  (in  the  sense 
that  a  small  change  in  entry  flight  path  angle  would  not  result  in  vehicle  loss),  and 
skip  trajectories  were  used  to  take  advantage  of  downrange  gains  incurred  by  such 
mission  profiles.  This  analysis  tended  toward  designs  with  base  cross  sections  that 
were  either  parallelograms  or  pure  ellipses  (or  combinations  of  the  two). 

As  before,  the  analytical  approach  under-predicts  DPLR  peak  stagnation  point 
heating  by  30%  to  70%,  and  the  aerodynamic  solutions  match  up  very  well  with  the 
lower-order  predictions  (within  10%  for  all  cases).  Cases  C  and  D  experience  their 
maximum  heat  flux  at  higher  altitudes  (above  64  km)  than  do  cases  A  and  F  (below 
60  km),  corresponding  to  both  the  lower-order  method  and  the  CFD  predicting  much 
lower  heating  rates  and  heat  loads.  Since  these  cases  do  their  primary  deceleration 
occurring  in  lower  density  atmosphere,  this  result  is  expected. 

The  following  subsections  detail  the  differences  between  the  the  lower-order 
methodology  and  DPLR  in  convective  heating  rates  for  the  different  cases  described 
in  Table  5.1.  Grid  topologies  are  not  shown  for  each  design,  but  all  CFD  meshes 
are  four-block  structured  grids  with  80  points  in  the  body  normal  direction.  Total 
grid  cells  for  each  case  are  tabulated  in  Table  5.2  at  the  end  of  this  section. 
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Table  5.1:  Summary  of  lunar  return  cases 


<3 

Q 


0) 

m 

cd 

o 


£ 

.2 

'C 

o 

<3 

o 


C3 

< 


S-i 

CD 

CD 

2 

Cd 

Jh 

cd 

Oh 


CO  ^  °  00  O 

CO  ;  _p  CO  O  H 
O  CM  CM  CM 


d4 


CO 


00 

o 

00 

CO 

o 

d4 

CO 

CO 

,co 

,co 

LO 

o 

CN 

'CN 

iO 

d4 

02 

o 

02 

o 

02 

o 

'CN 

iO 

02 

t-H 

'02 

3 

CO 

t-H 

'xf 

o 

LO 

00 

CN 

CO 

CO 

'cd 

ICO 

o 

t-H 

100 

,OI 

02 

t-H 

CD 

t-H 

cd 

CO 

CO 


o 

o 

C~^) 

o 

o 

00 

o 

lO 

,o 

00 

o  Ik. 

O  I  ■ 

3 

*o 

t-H 

CO 

o 

o 

t-H 

iO 

.•  '03 

<o 

CO 

'xf 

CN 

CN 

'xf 

o 

lO 

CN 

02 

t-H 

t-H 

CN 

'CO 

o 

t-H 

y-  Iio 

1  iCM 

iO 

t-H 

o 

t-H 

cd 

o  fr  o  CT3  03 

^  n  o  ro  b  1  1 

|m  n  H  o  12  s  S 


cm 

o 

Eh 

PL, 

Id 

'x 


db 


^0 


4-  Co 


kO 

Jh 

0) 

a 

o 

cd 

o 


t-h  CN  o 
lO  CM  00 

<9  ^  06 

O  O  CM 
LO  t-h  i 


CN  O 
CN  iO 


CD  Dt4  n. 


iO 


g 

’o 

Oh 

#bj) 

*c/2 

CD 

Q 


CO  02 


CN  of1  of 

b-  02 
LO  CO 


b-  b- 


o 

O  o 


00  lO  ^ 

o  h  c) 

T— I  T— I  (O’ 


e 

o 

<3 

J 


0^ 


o  qQ 


ptf 

HP 

PL, 

Q 


C/2 

CD 


cd 
Pi 
kO 

^  a 

CD  H 

-a  s 

C0  (X) 

3  ^ 

r— H  O 

Cd  H-l 


Dt4  db 

-d4 


01  >0> 


O  _ 

oo 


Dt4  b- 
O  T— I 

t-h  o 


o 
o 

d4  CN  m  CO 
^  lO  TO  ^ 

CM  t-h 


db  o  o  o 

ogro^ocoWo 

CM  - — -  CT3  O  O 

CM  CM  t-h I  t— I 

O  LQ  CO 

co 

H  O 


^  ^  ^ 

‘“d)  X-  oo  o  o 

00  <b3  CxO^-1  03  o  o 
CO  I  N  O  CO  o 

^  02  b~  T - 1  T - 1 

o  ^  CM 

^H  O 


c>  o 
—  o  o 


sis 

g  M  ig  g  x 

CO  1 — I  ° 


%  a 

S  CJ 
CD 

«  g 

0  5 


Q  Q 

^  s 


la  J 

«  e 
H  3 
e?  "s 
^  a. 


\  bO 
1  3  PP 


Q> 


c n 

o 

a 

cd 

~d  a 

(D  >2 

■+n  rx3 

o 


T3 

cd  h 

ft1  CD 

Ph  ^rH 

+3 

o 

$H 

CD 

c 


is 

II 

-I-H 

C/2 

c/2  0) 

CD 

CO 

c/2  C/2 

2 

GJ  rj 


SH 

?H  0) 

>  a 

g  ^ 

bJO  5P 

pi  3 
’w  s 

0  ? 

oT  qT 
be  W5 
Ed  a 
cd  cd 

Eh  '— h 

to  Ed 

TO  > 

o  § 

Eh  O 

o 

bO  bO 

.a  ^ 

N  ‘n 

a  a 

a’S 

a  a 

It3  03 
pd  pd 
cd  cd  yd 

N  ^  a 

Mg 

3  3 

<D 

-H  c^ 

Co  cd  _ 

1^  ^  8 

.a  a  I 

.2  'n  "a 

a  a  ” 


dd  pd 


CO 

CO 


85 


5.1.2  Case  A 


Case  A  (see  Figure  5.1)  is  a  slender  heat  shield  with  a  rounded  parallelogram 
base  much  like  those  discussed  in  Chapter  4.  A  shape  like  this  one  is  very  different 
from  a  classic  spherical  segment  and  may  be  more  difficult  to  implement  into  an 
actual  vehicle,  but  its  high  L/D  allows  it  to  possess  the  greatest  downrange,  or 
the  maximum  horizontal  distance  the  craft  travels  after  entry  interface,  of  all  cases. 
This  geometry’s  relatively  low  reference  area  ( Sref )  means  that  it  must  decelerate 
in  higher  density  atmosphere  thus  creating  the  most  adverse  heating  environment. 
Low  surface  area  corresponds  to  a  low  drag  area  ( Cc,Sref )  which  is  proportional  to 
drag  divided  by  dynamic  pressure,  itself  a  function  of  altitude  (free-stream  density). 
Basically,  to  achieve  the  same  amount  of  deceleration  (drag)  using  a  shape  with  a 
smaller  surface  area,  the  dynamic  pressure  must  be  higher.  This  is  achieved  only  at 
lower  altitudes  (corresponding  to  higher  free-stream  densities). 

DPLR  solutions  show  that  absolute  peak  heating  occurs  along  the  leading  edge 
of  the  vehicle  away  from  the  plane  of  symmetry  at  1420  W/cm2,  which  is  2.6  times 
the  stagnation  point  value  and  well  above  Orion  CEV  feasibility  limits.  Essentially, 
in  order  to  produce  a  design  with  greater  aerodynamic  maneuverability,  the  vehicle 
would  need  to  experience  heating  rates  higher  than  even  the  most  conservative 
estimates  for  Apollo.  The  situation  here  is  similar  to  what  was  observed  in  Chapter 
4.  The  lower-order  approach  seems  to  not  account  for  regions  of  potential  high  heat 
fluxes  away  from  the  stagnation  point.  These  results  further  emphasize  the  dangers 


in  using  a  parallelogram  base  cross  section  with  n2  close  to  1.3. 
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Figure  5.1:  Case  A  surface  pressure  and  convective  heating  profiles 
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5.1.3  Case  C 


Case  C  (see  Figure  5.2)  is  similar  in  shape  to  the  baseline  Orion  CEV  heat 
shield  (5.03  m  diameter,  9S  =  25.0°,  and  no  eccentricity).  This  geometry  is  shown 
here  to  provide  a  basis  of  comparison  for  the  other  optimized  designs.  This  design 
has  a  relatively  low  L/D ,  giving  it  the  lowest  cross  range  capabilities  of  all  designs 
studied  for  lunar  return.  Like  for  Apollo  4,  peak  convective  heating  occurs  on  the 
symmetry  plane  further  windward  of  the  stagnation  point  at  260  W/cm2.  The 
Orion  capsule’s  convective  wall  heat  flux  is  lower  than  that  of  Apollo  because  a 
larger  planform  area  allows  it  to  decelerate  much  higher  in  the  atmosphere  (64.1 
km  vs.  61  km).  Basically,  the  larger  drag  area  yields  a  lower  free  stream  density 
at  peak  instantaneous  heating,  and  the  lower  density  corresponds  to  lower  peak 
heating  rates. As  before,  the  maximum  convective  heating  pulse  is  1.66  times  the 
heating  experienced  at  the  stagnation  point  and  is  the  lowest  of  all  cases  examined 
here.  At  least  in  terms  of  convective  heating,  this  simple  geometry  would  appear 
to  provide  the  ideal  performance.  It  remains  to  be  seen,  however,  whether  or  not 
the  radiating  shock  layer,  necessary  to  produce  lower  convective  heating  rates,  will 
cause  total  heat  fluxes  to  exceed  what  is  experienced  by  the  other  designs. 
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Shape: 
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Figure  5.2:  Case  C  (Orion)  surface  pressure  and  convective  heating  profiles 


5.1.4  Case  D 


Case  D  (see  Figure  5.3)  is  an  oblate  design  with  an  L/D  slightly  larger  than 
that  of  Orion  («  0.3).  Here,  peak  heating  is  289  W/cm2,  which  is  1.49  times  the 
heating  at  the  stagnation  point.  Peaking  heat  is  spread  out  all  along  the  windward 
edge  of  the  heats  shield,  suggesting  that  the  effect  of  adding  edge  radius,  in  the 
form  of  a  torus,  is  to  temper  velocity  gradients  as  the  flow  is  turned  around  that 
edge.  Without  the  added  curvature,  peak  heating  rates  would  be  extremely  higher, 
leading  to  the  necessity  to  use  more  sophisticated  (and  more  expensive)  thermal 
protection  material.  This  case  also  exhibits  the  largest  relative  spread  between 
calculated  convective  heating  rates  and  low-order  predictions. 

The  lower-order  approach  suggested  that  this  geometry  would  experience  50% 
of  the  convective  heat  flux  experienced  by  Orion.  The  reason  being  that  this  design 
would  decelerate  at  a  higher  altitude  in  less  dense  atmosphere.  DPLR  results  show 
that  this  trend  does  not  actually  exist;  and  that,  in  fact,  this  design’s  stagnation 
point  heat  flux  actually  exceeds  that  calculated  for  the  Orion  analog.  This  obser¬ 
vation  supports  the  previous  assertion  that  non-axisymmetric  shapes  (i.e.  those 
generated  by  an  eccentric  base)  cause  adverse  heating  conditions  that,  in  turn,  force 
the  semi-empirical  correlations  used  in  the  lower-order  method  to  fail.  Without  flight 
data  for  such  eccentric  shapes,  it  is  nearly  impossible  to  discern  the  true  relationship 
between  eccentricity  and  the  resulting  heating  environment. 
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Figure  5.3:  Case  D  surface  pressure  and  convective  heating  profiles 


91 


5.1.5  Case  F 


Case  F  (see  Figure  5.4)  is  a  prolate  design  with  a  rounded  diamond  base 
optimized  using  a  conservative  mass  estimation  that  includes  a  three-fold  increase  in 
heat  shield  mass  to  deal  with  heating  loads.  Peak  heating  is  573  W/crn2,  1.45  times 
the  stagnation  point  rate;  and,  it  can  be  found  further  windward  on  the  symmetry 
plane  than  the  stagnation  point.  Greater  mass  forces  the  vehicle  to  decelerate  lower 
in  the  atmosphere,  yielding  high  convective  heating  rates,  but  not  more  than  the 
high  L/D  case. 

For  this  design,  the  low-order  stagnation  point  convective  heat  flux  under¬ 
predicts  the  DPLR  result  by  67%.  Both  this  case  and  the  previous  one  show  large 
discrepancies  in  predicting  the  stagnation  point  convective  heating  rates  using  the 
low-order  approach.  The  elliptical  nature  of  these  geometries  would  appear  to  force 
the  semi-empirical  correlations  to  report  incorrect  estimates.  Either  the  method 
in  which  effective  nose  radius  (the  driver  for  the  convective  heat  flux  relations)  is 
calculated  is  the  source  of  this  error  or  the  correlations  themselves  fail  to  account 
for  the  true  physical  nature  of  the  flow  around  such  shapes. 
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Figure  5.4: 


Case  F  surface  pressure  and  convective  heating  profiles 
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5.1.6  Computational  Cost  Summary 


Table  5.2  shows  a  summary  of  grid  sizes,  iteration  counts,  and  computational 
time  for  all  cases  considered  in  this  section.  No  case  took  less  than  30,000  iterations 
and  1,400  CPU  hours  to  converge.  Case  A  needed  the  most  iterations  to  converge 
to  a  stable  solution  due  to  a  small  pocket  of  subsonic  flow  existing  at  the  windward 
exit  of  its  mesh.  More  care  was  needed  to  ensure  that  this  case  was  not  influenced  by 
this  discrepancy  and,  in  fact,  did  reach  a  stable  solution.  There  exists  a  possibility 
that  this  boundary  condition  violation  would  introduce  errors  into  the  final  solution 
for  this  case;  however,  results  for  that  design  are  consistent  with  a  similar  shape 
(see  Chapter  4),  suggesting  that  the  errors,  if  they  exist,  are  negligible. 


Table  5.2:  Summary  of  computational  costs  for  lunar  return  cases 


Case 

Grid  Size  (#  cells) 

Iterations 

CPU  Time  (hrs) 

A 

268800 

85400 

3247.16 

C 

250880 

47100 

1882.16 

D 

336000 

49300 

3640 

F 

232960 

30000 

1442.78 

94 


5.2  Mars  Return  Optimized  Designs 


This  section  presents  computational  solutions  for  heat  shields  generated  for 
earth  return  after  a  mission  to  Mars  using  the  coupled  trajectory/ vehicle  optimiza¬ 
tion  scheme  discussed  in  Chapter  1.  Again,  no  NEQAIR  simulations  are  included 
for  these  cases  due  to  time  and  computational  constraints.  As  such,  only  convective 
heating  and  aerodynamic  predictions  are  presented.  A  possible  breakdown  of  the 
continuum  flow  assumption  used  by  DPLR  is  discussed  in  this  section,  but  non¬ 
continuum  simulations  are  left  for  future  work. 

5.2.1  General  Summary 

DPLR  results  are  obtained  for  heat  shields  generated  using  the  coupled  op¬ 
timization  technique  for  Mars  return  entry  velocities,  Ve  =  12.5  krn/s,  and  entry 
flight  path  angles  of  =  —6.4°  at  the  location  on  the  trajectory  upon  which  the 
predicted  peak  instantaneous  heat  flux  occurs.  Tabic  5.3  shows  a  summary  of  the 
geometry,  design  point,  the  aerothermodynamics  calculated  using  DPLR,  and  the 
lower-order  predicted  aerothermodynamics  for  all  cases  studied  in  this  section.  The 
cases  maintain  their  descriptors  from  Table  13.1  of  Reference  [6].  Percent  differ¬ 
ences,  in  reference  to  the  DPLR  solutions,  are  presented  for  the  lower-order  results 
in  parenthesis  where  applicable.  Both  spherical  segment  and  sphere-cone  axial  pro¬ 
files  are  considered  here  as  the  optimizer  generated  independent  geometries  for  these 
two  types  of  topologies.  All  power-law  optimized  shapes  were  simply  disguised  ver¬ 
sions  of  the  two  other  profiles.  Optimizations  were  performed  at  L/ D  =  0.3  and  0.5 
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for  two  objective  functions  sets:  (1)  maximizing  downrange  while  minimizing  heat 
load  and  (2)  maximizing  cross  range  while  minimizing  heat  load.  No  optimizations 
were  done  at  L/D  =  1.0  because  those  designs  resulted  in  heat  loads  that  were 
considered  infeasible.4  Smaller  entry  corridor  widths  of  up  to  0.79°  were  necessary 
to  generate  flyable  trajectories  here. 

Similarly  to  the  lunar  return  cases,  designs  with  the  lowest  heat  loads  will 
experience  their  peak  heat  pulses  at  much  higher  altitudes.  For  all  Mars  return 
cases,  DPLR  reports  stagnation  point  convective  heat  fluxes  much  lower  than  does 
the  simple,  analytical  method.  For  example,  the  stagnation  point  heat  flux  for  case 
B  is  only  19%  that  of  what  was  calculated  using  the  modified  Newtonian  approach. 
At  Mars  return  velocities,  the  shock  layer  is  actually  larger  than  at  lower  speeds 
(larger  shock  stand  off  distance);  and,  consequently,  shock  layer  radiation  should 
play  a  larger  role  in  the  resulting  heating  profile  for  a  given  blunt-body  heat  shield. 
A  firm  grasp  of  the  potentially  strong  coupling  between  convection  and  radiation 
must  be  reached  before  making  any  concrete  conclusions  about  the  designs  studied 
in  this  section.  Still,  there  is  a  great  deal  to  glean  from  comparing  the  higher 
order  simulation  with  lower-order  predictions,  especially  if  the  aim  is  to  improve  the 
lower-order  method  for  use  in  initial  design  studies. 

Furthermore,  most  of  these  cases  require  the  addition  of  increased  numerical 
dissipation  for  convergence.  This  increased  dissipation  is  necessary  due  to  the  possi¬ 
bility  of  non-continuum  flow  present  at  these  entry  conditions.  Continuum  flow  can 
be  classified  through  the  use  of  the  Knudsen  number,  Km,  which  is  the  ratio  of  mean 
free  path  (distance  a  molecule  will  travel  before  colliding  with  another  molecule)  to 
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a  characteristic  length.  Mean  free  path,  in  turn,  is  a  function  of  temperature  divided 
by  pressure.13  Regions  of  relatively  high  temperature  and  relatively  low  pressures 
are  more  likely  for  the  higher  Mars  return  velocities,  yielding  values  of  that 
may  violate  the  continuum  requirement  of  KN  <  0.3.  Another  way  to  classify 
non-continuum  flow  is  by  using  the  gradient- length  local  Knudsen  number:73 


K 


N,GLL 


X 

Q 


dQ 


dl 


(5.1) 


where  A  is  the  mean  free  path,  Q  is  a  flow  property  (usually  temperature),  and  l  is  a 
distance  between  two  points  in  the  flow  field  along  the  direction  of  steepest  gradients. 
Continuum  breakdown  occurs  when  the  value  of  this  parameter  is  less  than  0.05, 
and  this  definition  of  Knudsen  number  is  better  suited  to  the  discretized  flow  holds 
used  in  computational  fluid  dynamics  as  it  is  relatively  easy  to  extract  gradients 
from  a  computational  solution.  Adding  extra  dissipation  may  force  continuum  how 
to  exist  everywhere  in  the  how  held,  but  it  adds  further  sources  of  error  to  the 
solutions.  Though  results  are  consistent  with  what  is  seen  at  lunar  return  velocities, 
it  is  paramount  that  this  additional  source  of  error  be  classihed  and  quantified  before 
robust  conclusions  are  made.  Such  classifications  are  left  in  the  realm  of  future  work. 

The  following  subsections  detail  the  differences  between  the  the  lower-order 
methodology  and  DPLR  in  convective  heating  rates,  as  the  aerodynamic  predictions 
are  almost  identical,  for  the  different  cases  described  in  Table  5.3.  Grid  topologies 
are  not  shown  for  each  design,  but  all  CFD  meshes  are  four-block  structured  grids 
with  80  points  in  the  body  normal  direction.  Total  grid  cells  for  each  case  are 
tabulated  in  Table  5.4  at  the  end  of  this  section. 
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Table  5.3:  Summary  of  Mars  return  cases 
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5.2.2  Case  A 


Case  A  (see  Figure  5.5)  is  a  prolate  shape,  similar  to  Case  F  for  lunar  return, 
with  the  highest  L/D  of  all  shapes  considered  here;  however,  this  maneuverability 
comes  at  the  cost  of  having,  by  far,  the  highest  maximum  convective  heat  flux.  Peak 
convective  heating  (657  W/cm2)  occurs  along  the  axis  symmetry  further  windward 
of  the  stagnation  point  and  is  1.47  times  greater  than  the  calculated  value  there. 
The  lower-order  stagnation  point  convective  heat  flux  is  64%  lower  than  the  DPLR 
calculation,  which  is  the  largest  discrepancy  in  this  category  for  all  cases  for  Mars 
return. 

Notably,  this  case  is  the  only  one  studied  in  this  section  that  does  not  need 
extra  numerical  dissipation  to  converge.  This  suggests  that  continuum  flow  assump¬ 
tion  is  valid  for  this  design  at  its  predicted  peak  heating  trajectory  point  and  that 
the  conclusions  drawn  from  these  results  are  free  of  the  additional  errors  that  plague 
the  other  heat  shields  examined  for  Mars  return.  Case  F  for  lunar  return,  which  is 
nearly  the  same  shape  as  this  design,  shows  an  almost  identical  offset  for  stagnation 
point  convective  heat  flux  (a  67%  difference  between  the  lower-order  methodology 
and  DPLR).  In  that  respect,  this  severe  under-prediction  by  the  lower-order  ap¬ 
proach  implies  that  the  semi-empirical  correlations  are  not  suited  to  a  shape  of  this 
class  for  the  same  reasons  as  discussed  in  Section  5.1.5. 
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Figure  5.5:  Case  A  surface  pressure  and  convective  heating  profiles 
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5.2.3  Cases  B  and  F  -  Orion  Analogs 


Cases  B  (see  Figure  5.6)  and  F  (see  Figure  5.7)  are  Orion  sized  spherical 
segments  (5.03  m  diameter,  9S  =  25.0°,  no  eccentricity)  simulated  using  a  lower 
(no  additional  heat  shield  mass  added  to  account  for  heat  loads)  and  an  upper 
(three-fold  increase  in  heat  shield  mass)  mass  estimation  respectively.  Both  cases 
have  peak  heating  occurring  along  the  axis  of  symmetry  further  windward  of  the 
stagnation  point.  Peak  heat  flux  is  1.45  (320  W/cnr)  and  1.75  (361  W/cm2)  times 
that  experienced  at  the  stagnation  point  for  the  lower  and  upper  mass  estimations 
respectively.  These  particular  values  bracket  the  baseline  Apollo  value  (1.66)  and 
fall  well  within  the  span  of  what  was  observed  for  lunar  return. 

Both  modified  Newtonian  solution  sets  report  stagnation  point  convective  heat 
fluxes  that  are  relatively  similar  to  their  DPLR  counterparts  (7%  and  23%  less 
respectively).  Since  added  dissipation  was  needed  for  these  cases,  it  is  reasonable  to 
assume  that  the  true  difference  between  the  two  approaches  will  actual  be  larger  (on 
the  order  of  what  was  observed  for  Apollo  and  axisymmetric  lunar  return  designs). 
When  mass  is  increased,  the  lower-order  method  shows  a  stagnation  point  heat 
flux  increase  of  12%.  Essentially,  additional  mass  translates  to  a  lower  altitude 
deceleration  and,  thus,  higher  heating  rates.  DPLR  solutions,  on  the  other  hand, 
show  a  36%  increase  convective  heating  rate  at  the  stagnation  point.  This  is  a 
relatively  modest  increases  that  may  be  an  artifact  pf  the  added  dissipation.  Still, 
nothing  in  the  CFD  results  suggest  that  anything  other  than  altitude  is  responsible 
for  the  increase  in  convective  heat  flux. 


101 


Figure  5.6:  Case  B  (Orion)  surface  pressure  and  convective  heating  profiles 
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Figure  5.7:  Case  F  (Orion)  surface  pressure  and  convective  heating  profiles 
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5.2.4  Case  D 


Case  D  (see  Figure  5.8)  represents  an  oblate  design,  with  similar  L/D  to  the 
Orion  analog  heat  shield,  optimized  for  minimized  heat  load  an  maximum  cross 
range  using  the  lower  mass  estimation.  Peak  heating  is  spread  out  along  the  wind¬ 
ward  edge  of  the  vehicle  upstream  of  the  stagnation  point.  Peak  convective  heat 
flux  (236  W/cm2)  is  2.15  times  more  that  which  is  experienced  at  the  stagnation 
point  (128  W/cm2),  and  the  lower-order  approach  under-estimates  stagnation  point 
heating  rates  by  approximately  30%.  Though  this  case  has  an  eccentric  base,  the 
discrepancy  between  DPLR  and  low-order  results  is  on  the  order  of  what  was  ob¬ 
served  for  Apollo  4.  Simply,  while  eccentricity  pushes  the  convective  heat  flux  up, 
the  added  dissipation  drops  the  calculated  rate,  creating  a  false  sense  of  consistency 
with  axisymmetric  cases  (the  lower  order  method  under-predicts  computational  so¬ 
lutions  by  ~  70%  for  all  other  elliptical  base  cases). 

This  case  was  meant  to  represent  a  marked  improvement  over  the  baseline 
Orion  geometry  as  its  larger  surface  area  should  allow  for  lower  heat  loads  and  heat 
fluxes,  while  maintaining  similar  aerodynamic  performance  (i.e.  Lj  D ,  pxrs,  Pdwn )• 
At  first  glance,  the  DPLR  results  show  such  an  improvement.  A  comparison  of  the 
two  cases  reveals  that  peak  convective  heating  decreases  by  a  modest  26%  while 
stagnation  point  convective  heating  lowers  by  30%  when  comparing  this  case  to 
the  low  mass  Orion  analog.  Similarly,  the  modified  Newtonian  approach  predicts  a 
slightly  larger  35%  decrease  in  stagnation  point  heat  flux  between  the  two  cases.  The 
more  elliptical  heat  shield  shows  a  marked  decrease  in  peak  heating  when  compared 
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to  a  simple  spherical  segment.  This  is  the  opposite  trend  than  was  observed  for 
lunar  return.  There,  the  heating  rate  actually  increased  when  transitioning  from 
an  axisymmetric  heat  shield  to  an  elliptical  one.  The  disparity  that  exists  between 
the  different  flight  regimes  stems,  most  likely,  from  the  poor  classification  of  non¬ 
continuum  effects  that  precipitates  the  need  for  added  numerical  dissipation,  rather 
than  something  physically  different  in  the  flow  fields  for  Mars  and  lunar  return. 

Furthermore,  it  is  interesting  to  note  that  surface  heating  contours  for  Case  D 
show  that  peak  heating  is  more  spread  out  over  the  entire  windward  edge  of  the  heat 
shield  as  opposed  to  the  more  local  and  concentrated  heat  pulse  displayed  in  case 
B.  This  means  that  the  highest  heat  loads  would  potentially  be  more  widely  spread 
over  the  elliptic  heat  shield,  forcing  the  addition  of  more  thermal  material  which,  in 
turn,  adds  to  vehicle  weight.  Simply,  there  is  no  single  metric  here  that  can  prove 
whether  or  not  this  shape  is  really  an  improvement  over  the  simpler  axisymmetric 
geometry.  Really,  until  all  the  errors  of  non-continuum  flow  and  the  elliptical  effects 
can  be  quantified  into  improved  empirical  correlations,  it  will  always  be  difficult  to 
determine  which  shape  is  more  ideal. 
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Figure  5.8:  Case  D  surface  pressure  and  convective  heating  profiles 
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5.2.5  Computational  Cost  Summary 


Table  5.4  shows  a  summary  of  grid  sizes,  iteration  counts,  and  computational 
time  for  all  cases  considered  in  this  section.  No  case  took  less  than  30,000  iterations 
and  800  CPU  hours  to  converge.  Case  D  needed  the  most  time  to  arrive  at  a  steady 
state  solution,  possibly  due  to  side-effects  from  the  added  dissipation.  Case  A,  which 
needed  no  extra  numerical  dissipation,  converged  the  fastest,  while  the  other  cases 
were  all  much  more  computationally  intensive.  Probably,  the  numerical  dissipation, 
though  helping  keep  the  solution  stable,  is  the  source  of  the  slow  convergence  rates. 


Table  5.4:  Summary  of  computational  costs  for  Mars  return  cases 


Case 

Grid  Size  (#  cells) 

Iterations 

CPU  Time  (hrs) 

A 

219520 

30500 

847.84 

B 

250880 

52600 

1773.38 

D 

250780 

91100 

2714.46 

F 

250880 

40600 

1446.67 
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Chapter  6 


Conclusions 

6.1  Summary  of  Results 

In  general,  computational  fluid  dynamics  solutions  of  the  Apollo  4  capsule  and 
of  optimized  heat  shield  geometries  show  that  the  lower-order  approach  discussed  in 
Chapter  1  gives  reasonable  estimations  of  aerothermodynamic  properties  useful  for 
initial  design  studies.  Particularly,  CFD  solutions  show  that  the  modified  Newtonian 
approach,  as  expected,  gives  highly  accurate  predictions  for  the  aerothermodynamic 
parameters  (i.e.  Cl,  Cd,  and  L/D );  however,  large  disparities  in  convective  and 
radiative  (where  applicable)  heating  profiles  are  seen.  The  following  subsections 
detail  the  important  results  discussed  in  this  work. 

6.1.1  Apollo  4  Benchmarking 

For  the  Apollo  axisymmetric  heat  shield  at  Apollo  4  peak  heating  conditions, 
the  lower-order  approach  under-predicts  convective  stagnation  point  flux  by  30% 
and  over  predicts  stagnation  point  radiative  heat  flux  by  16%  when  compared  to 
computational  solutions.  These  disparities  can  be  attributed  to  the  failure  of  the 
lower-order  method  to  capture  boundary  layer  physics.  The  correlations  used  to 
calculate  the  heating  profiles  in  the  analytical  approach  were  formulated  for  simple 
axisymmetric  shapes  like  spheres;  as  such,  these  offsets  provide  a  useful  baseline 
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for  all  other  comparisons  of  the  computational  solutions  to  the  lower-order  ones. 
Essentially,  if  the  errors  seen  for  other  geometries  exceed  the  ones  seen  here,  then 
other  inaccuracies,  above  and  beyond  boundary  phenomena,  must  exist  in  the  lower- 
order  approach.  Though  these  errors  are  modest  but  significant,  combining  the  two 
heating  regimes  yields  an  offset  of  only  8%  for  total  heat  flux  between  the  two 
approaches.  The  errors  that  manifest  themselves  in  the  low-order  correlations  used 
to  calculate  heat  transfer  would  appear  to  cancel  each  other  out. 

Edge  radius  plays  a  major  role  in  the  aerothermodynamics  of  a  blunt-body 
heat  shield,  especially  with  regards  to  its  heating  environment.  Ignoring  curvature 
at  the  edges  of  the  heat  shield  introduces  discontinuities  in  the  heating  profile  and 
may  even  cause  a  violation  of  DPLR  boundary  conditions.  For  the  Apollo  capsule, 
peak  heating  decreases  as  a  power  law  function  of  the  exact  curvature  that  exists 
at  the  edge  of  the  vehicle,  even  at  an  angle  of  attack  that  would  place  the  stagna¬ 
tion  point  further  away  from  the  windward  edge  of  the  heat  shield.  Including  edge 
radius  as  a  design  variable  in  the  design  process  may  prove  difficult  without  more 
detailed  mission  profiles,  but  is  necessary  to  add  this  feature  if  the  a  true  optimum 
geometry  is  desired.  One  way  to  implement  edge  curvature  into  the  optimization 
process  would  be  to  blend  a  torus,  of  either  fixed  or  variable  radius,  to  the  geome¬ 
tries  generated  by  the  process  discussed  in  Chapter  1.  This  addition  would  cause 
additional  computational  cost,  through  the  addition  of  mesh  points  and  the  possible 
inclusion  of  more  optimization  constraints,  but  the  results  would  provide  a  much 
more  accurate  representation  of  what  would  be  expected  aerothcrmodynamically 
from  an  actual  blunt-body  heat  shield. 
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6.1.2  High  L/D  shapes 

For  a  representative  slender,  high  L/D  heat  shield  with  n2  =  1.3,  the  lower- 
order  approach  under-estimates  stagnation  point  convective  heating  by  44%  and 
over-predicts  stagnation  point  radiative  heat  flux  by  79%  compared  to  the  compu¬ 
tational  solutions.  For  convective  heat  flux  the  percent  offset  is  similar  to  what  was 
observed  for  the  Apollo  benchmark  case  because  the  stagnation  point  falls  in  a  highly 
spherical  region  of  the  heat  shield.  Still  the  14%  increase  (from  30%  seen  for  Apollo 
4  to  44%  seen  here)  in  convective  heat  flux  offset  (comparing  the  lower-approach  to 
DPLR  solutions)  suggests  that  the  elliptical  nature  of  the  base  cross  section  may 
cause  the  relations  sued  to  predict  this  value  in  the  lower-order  approach  to  break 
down.  The  stagnation  point  radiative  heat  flux  offset  is  much  larger  than  the  16% 
seen  for  the  axisymmetric  Apollo  case.  Likely,  the  process  by  which  shock  stand-off 
distance  (the  driving  factor  for  radiative  heat  flux)  is  calculated  is  not  suitable  for 
a  slender  body  such  as  this  one.  Still,  more  evidence  to  this  effect  would  need  to  be 
accrued  before  this  assertion  could  be  truly  corroborated.  Again,  combining  the  two 
heating  rates  calculated  computationally  generates  a  result  that  compares  favorably 
(within  20%)  to  what  was  predicted  using  the  lower-order  approach. 

At  Apollo  4  peak  heating  conditions,  using  ri2  =  1.3  does  indeed  produce 
high  lift  geometries;  however,  the  heating  profiles  for  these  shapes  show  maximum 
convective  heating  to  be  more  than  twice  what  is  experienced  at  the  stagnation 
point.  The  effective  sharpness  of  the  base  cross-section  creates  a  sharp  leading  edge 
near  the  boundaries  of  the  heat  shield  when  the  axial  profile  is  added  to  complete 
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the  geometry.  While  not  exactly  a  discontinuity,  the  low  bluntness  at  this  leading 
edge  generates  a  thin  shock  layer  that  contributes  to  the  adversely  high  heating  rates 
experienced  near  the  edges  of  the  heat  shield.  This  same  thin  shock  layer  generates 
negligible  radiative  heat  transfer  at  the  point  of  highest  heating,  but  the  maximum 
convective  heating  is  still  almost  twice  that  of  the  total  heat  flux  at  the  stagnation 
point  (where  peak  radiative  heating  occurs).  Some  newer  materials  might  be  able 
to  withstand  heating  rates  at  or  near  1000  W/cm2,  but  that  would  push  design 
limitations  imposed  for  the  current  Orion  CEV  capsule  and  significantly  add  to 
vehicle  cost.  For  the  representative  high  L/D  shape,  maximum  convective  heating 
decreases  as  a  power  law  as  the  712  parameter  is  increased.  A  an  approximate  40% 
reduction  in  peak  convective  heating  can  be  obtained  by  increasing  the  ri2  parameter 
to  1.5  while  still  maintaining  an  L/D  greater  than  1.  As  such,  future  optimizations 
should  alter  the  lower  bound  on  this  n,2  parameter  to  1.4  or  1.5  in  order  to  generate 
high  L/D  geometries  without  the  adverse  off  stagnation  point  heating  seen  here.. 

6.1.3  Coupled  Vehicle/Trajectory  Optimized  Geometries 

For  lunar  return,  shapes  with  eccentric  bases  (either  prolate  or  oblate)  show 
qualitative  discrepancies  in  heating  profile  when  compared  to  the  modified  New¬ 
tonian  solutions.  A  close  examination  of  DPLR  solutions  show  that  any  possible 
gains  from  increased  surface  area  and  higher  altitude  decelerations  are  wiped  out  by 
changes  in  the  flow-field  introduced  by  stretched  geometries.  As  such,  there  is  reason 
to  suspect  that  these  eccentric  shapes  fall  outside  the  realm  of  the  semi-empirical 
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correlations  used  in  the  lower-order  method  for  heat  fluxes.  Simply,  the  correlations 
are  not  valid  for  every  class  of  blunt-body  heat  shield  created  by  the  lower-order 
optimization  process.  These  relations  can,  however,  be  improved  through  the  use 
of  a  larger  data  set  of  wind  tunnel  data  and  CFD  solutions  (like  those  seen  in  this 
thesis)  that  includes  more  shapes  with  eccentric  bases  and  sharp  edges.  For  Mars 
return,  errors  associated  with  non-continuum  flow  (manifesting  itself  in  increased 
dissipation)  and  eccentricity  effects  make  it  difficult  to  make  any  concrete  conclu¬ 
sions  regarding  the  merits  of  one  design  over  another.  These  errors  need  to  be 
quantified  and  accounted  for  before  any  such  study  may  continue. 

Furthermore,  it  can  be  seen  that  it  is  difficult  to  produce  heat  shields  that 
show  a  great  deal  of  improvement,  for  both  lunar  and  Mars  return,  over  one  with 
geometric  parameters  similar  to  that  which  is  currently  in  consideration  for  the 
Orion  CEV  capsule.  Any  advantages  gained  by  using  a  more  novel  shape  will,  more 
than  likely,  be  canceled  out  by  the  ease  of  manufacturing  and  vehicle  integration  for 
an  Apollo-like  spherical  segment  design.  Whether  intended  or  not,  it  would  appear 
that  a  simple  spherical  segment  with  6S  =  25.0°  is  indeed  an  ideal  shape  for  earth 
entry  at  super  orbital  velocities. 

6.2  Future  Work 

Future  additions  to  this  work  fall  into  four  categories:  1)  better  radiation 
modeling,  2)  material  response,  3)turbulcnce,  and  4)  other  atmospheres  for  entry. 
The  heating  profile  for  a  blunt-body  heat  shield  can  not  accurately  be  calculated 
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unless  convection  and  radiation  are  fully  coupled.  Loosely  coupling  the  two  with 
NEQAIR  and  DPLR  is  possible  with  the  technique  described  in  Chapter  2,  but 
future  versions  of  DPLR  will  contain  internal  methods  for  dealing  with  shock  layer 
radiation,  allowing  for  much  simpler  acquisition  of  fnlly  coupled  solutions  without 
the  high  computation  costs  associated  with  using  NEQAIR.  Materials  play  a  ma¬ 
jor  role,  through  ablation,  in  determining  the  heat  actually  felt  by  the  vehicle.  All 
materials  will  undergo  sometime  sort  of  chemical  change  when  exposed  to  the  ex¬ 
treme  environments  experienced  during  re-entry.  Future  CFD  solutions,  and  the 
lower-order  optimized  geometries  for  that  matter,  need  to  take  into  account  how 
chemical  changing  in  a  vehicle’s  surface  will  change  the  resulting  flow-field  around 
a  next  generation  space  capsule  (i.e.  the  gas  model  changes  as  the  environment  is 
no  longer  just  air)  if  truly  accurate  solutions  are  desired. 

All  CFD  solutions  in  this  work  assume  a  laminar  flow.  It  is  not  entirely  obvious 
weather  or  not  earth  entering  heat  shields  will  experience  local  regions  of  turbulent 
flow.  To  that  end,  it  would  be  germane  to  adopt  a  some  sort  of  transition  criteria, 
based  on  flow  physics,  and  then  apply  turbulence  models,  within  DPLR,  to  those 
regions  to  correctly  model  the  flow.  This  transition  criteria  might  be  more  pertinent 
in  different  planetary  atmospheres  such  as  Mars,  where  turbulent  flow  is  much  more 
likely  to  exist.  Furthermore,  since  a  next  generation  space  capsule  will  be  used 
for  missions  that  require  entry  to  the  atmospheres  of  other  planets,  it  would  be 
interesting  to  see  how  the  shapes  studied  in  this  work,  optimized  for  earth  entry, 
measure  up  in  different  environments,  like  that  on  Mars. 
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6.3  Concluding  Remarks 


At  this  point,  it  is  important  to  try  to  understand  the  computational  results 
compiled  in  this  work  in  a  more  general  perspective.  The  purpose  of  conducting 
high  fidelity  CFD  on  geometries  optimized  using  a  low-order  approach  was  to  evalu¬ 
ate  how  well  that  analytical  method  could  predict  the  extreme  aerothermodynamic 
environments  these  geometries  would  experience  on  a  real  mission.  To  that  end, 
the  CFD  results  show  that  the  aerodynamic  model  used  by  the  low-order  approach 
does  a  more  than  adequate  job  in  predicting  the  proper  pressure  distribution  on 
the  heat  shield  surface,  while  the  correlations  used  to  predict  the  thermodynamic 
environment  prove  poor,  even  in  circumstances  for  which  they  were  derived  for  (i.e. 
spheres).  Also,  these  stagnation  point  heating  models  fail  to  pick  up  areas  of  high 
heat  flux  on  other  parts  of  the  vehicle’s  surface,  highlighting  a  further  shortcom¬ 
ing  of  the  lower-order  analytical  approach.  All  of  these  conclusions  were,  to  some 
extent,  expected;  but  the  process  by  which  they  were  obtained  implies  possible  im¬ 
provements  for  the  lower-order  approach.  Namely,  that  the  empirical  relations  used 
to  predict  heating  rates  need  to  be  replaced  with  improved  correlations  with  a  more 
physical  basis  and  that  the  geometric  constraints  used  in  the  optimization  process 
need  to  be  further  limited  in  order  to  avoid  large  local  off-stagnation  point  heat 
fluxes. 

Furthermore,  the  results  gathered  in  this  work  allow  for  some  conclusions 
to  be  made  about  the  blunt-body  design  space  in  general.  Presently,  there  are 
very  few  physical  data  points  for  blunt-body  entry  at  extra-planetary  velocities. 
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As  such,  the  empirical  relations  derived  from  those  data  points  break  down  when 
they  are  used  for  unconventional  shapes.  In  order  to  improve  these  correlations, 
more  data  must  be  obtained  through  flight  and  ground  tests  as  well  as  through 
further  computational  simulations.  That  way,  curve  fits  derived  from  this  larger 
data  set  will  truly  reflect  the  full  range  of  possible  outcomes.  Finally,  the  lower- 
order  method  sought  to  show  that  more  complicated  shapes  could  provide  large  gains 
over  the  simpler,  axisymmetric  geometries.  However,  in  practice,  the  more  complex 
shapes  introduce  aspects  into  the  blunt-body  flow  field  (i.e.  high  off  stagnation  point 
heating  and  other  elliptical  effects)  that  do  not  manifest  themselves  with  the  simpler 
shapes.  Sometimes  the  simpler  approach  can  be  the  better  one;  and,  certainly  in 
this  work,  it  can  be  seen  that  choosing  a  simpler  shape  (in  this  case  a  25°  spherical 
segment)  can  be  more  advantageous,  in  many  respects,  than  a  more  complicated 
design. 
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